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ABSTRACT 

We use hydrodynamical simulations of disk galaxies to study relations between star formation and properties 
of the molecular interstellar medium (ISM). We implement a model for the ISM that includes low-temperature 
(T < 10 4 K) cooling, directly ties the star formation rate to the molecular gas density, and accounts for the 
destruction of H2 by an interstellar radiation field from young stars. We demonstrate that the ISM and star for- 
mation model simultaneously produces a spatially-resolved molecular-gas surface density Schmidt-Kennicutt 
relation of the form Ssfr oc E^™ 1 with n mo \ w 1.4 independent of galaxy mass, and a total gas surface density - 
star formation rate relation Ssfr oc EJg with a power-law index that steepens from n t ot ~ 2 for large galaxies 
t° "tot ^ 4 for s mall dwarf galax ies. We show that deviations from the disk-averaged Ssfr oc Ei^ correlation 
determined by lKennicuttl (119981) owe primarily to spatial trends in the molecular fraction fm and may explain 
observed deviations from the global Schmidt-Kennicutt relation. In our model, such deviations occur in re- 
gions of the ISM where the fraction of gas mass in molecular form is declining or significantly less than unity. 
Long gas consumption time scales in low-mass and low surface brightness galaxies may owe to their small 
fractions of molecular gas rather than mediation by strong supernovae-driven winds. Our simulations also re- 
produce the observed relations between ISM pressure and molecular fraction and between star formation rate, 
gas surface density, and disk angular frequency. We show that the Toomre criterion that accounts for both 
gas and stellar densities correctly predicts the onset of star formation in our simulated disks. We examine the 
density and temperature distributions of the ISM in simulated galaxies and show that the density probability 
distribution function (PDF) generally exhibits a complicated structure with multiple peaks corresponding to 
different temperature phases of the gas. The overall density PDF can be well-modeled as a sum of lognormal 
PDFs corresponding to individual, approximately isothermal phases. We also present a simple method to miti- 
gate numerical Jeans fragmentation of dense, cold gas in Smoothed Particle Hydrodynamics codes through the 
adoption of a density-dependent pressure floor. 

Subject headings: Galaxies, Star Formation 



1. INTRODUCTION 

Galaxy formation presents some of the most important 
and challenging problems in modern astrophysics. A ba- 
sic paradigm for the dissipational formation of galaxies from 
primordial flu ctuations in the density field has been de- 
veloped (e.g. IWhite & Reesl 119781; iBlumenthal etail 11984 
IWhite & Frenklll991l) . but many of the processes accompa- 
nying galaxy formation are still poorly understood. In partic- 
ular, star formation shapes the observable properties of galax- 
ies but involves a variety of complicated dynamical, thermal, 
radiati ve, and chemical proce sses on a wide range of scales 
(see lMcKee & O striker 20071 for a review). Observed galax- 
ies exhibit large-scale correlations between their global star 
formation rate (SFR ) surface density Esfr and average gas 
surface density S gas dKennicuttll 1 9891 .1 1 9981) . and these global 
correlations serve as the basis for treatments of star forma- 
tion in many models of galaxy formation. While such mod- 
els have supplied important insights, detailed observations of 
galaxies have recently provided evidence that the molecular, 
rather than the total, gas surface density is the primary driver 
of global star format i on in galaxies (e.g . JWong & Blit j|2002 ; 



Boiss ier et all l2003t iHever et all l2004 iBoissier et al l 12007 
Calzetti et al.ll2007HKennicutt et al.ll2007l) . 
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In this study, we adopt an approach in which empirical and 
theoretical knowledge of the star formation efficiency (SFE) 
in dense, molecular gas is used as the basis for a star forma- 
tion model in hydrodynamical simulations of disk galaxy evo- 
lution. This approach requires modeling processes that shape 
properties of the dense phase of the interstellar medium (ISM) 
in galaxies. The purpose of this paper is to present such a 
model. 

Stellar populations in galaxies exhibit salient trends 
of colors and metallic i ties with galaxy luminosity (e.g., 
Kauffmannetai] 120031: iBlanton et all 120051: ICooper et al.1 
l2007t> . In the hierarchical structure formation scenario these 
trends should emerge through the processes of star forma- 
tion and/or stellar feedback in the progenitors of present-day 
galaxies. Observationally, ample evidence suggests that the 
efficiency of the conversion of gas into stars depends strongly 
and non-monotonically on mass of the system. For example, 
the faint-end of the galaxy lumin osity function has a shal- 
low slope (ol « 1.0-1.3, e.g.. IBlanton et al.1 1200 ll l2003h 
compared to the steeper mass function of dark matter ha- 
los (a DM « 2, e.g. JPress & Schechterlll974ISheth & Tormenl 
1999), indicating a decrease in SFE in low-mass galaxies. At 
the same time, the neutral hydrogen (HI) and baryonic mass 
functions may b e steeper than the luminosity fu n ction (am sa 
1.3 — 1.5, e.g., iRosenberg & Schneider! l2002t IZwaan et al.l 
120031) . The baryonic iTully & Fisherl (1 19771) relation is con- 
tinuous down to extremely low-m ass dwarf galaxies (e.g., 
iMcGaughl 120051: iGeha et al J 120061) . indicating that the frac- 
tional baryonic content of galaxies of different mass is similar. 



2 



Robertson and Kravtsov 



Hence, low-mass galaxies that are unaffected by environmen- 
tal processes are gas-rich, yet often form stars inefficiently. 

While feedback processes from super novae and AGN (e.g., 
iBrooks et alj 120071 : ISiiacki et al.l |2007[). or the efficiency o f 
gas cooling and accretion (jDekel & Birnboiml 120061 2008), 
may account for part of these trends, the SFE as a function 
of galaxy mass may also owe to intrinsic ISM processes (e.g., 
Tassis et al. 2008; Kauf mann et al.ll2~007l) . To adequately ex- 
plore the latter possibility, a realistic model for the conversion 
of gas into stars in galaxies is needed. 

Traditionally, star formation in numerical simulations of 
galaxy formatio n is based on the empirical Schmid t-Kennicutt 
(SK) relation (ISchmidtl 1 959HKennicuttl 1 98l[l998h . in which 
star formation rate is a universal power-law function of the to- 
tal disk-averaged or global gas surface density: Esfr oc Eg™ 
with n tot K 1.4 describing the correlation for the entire pop- 
ulation of normal and starburst galaxies. However, growing 
observational evidence indicates that this relation may not be 
universal on smaller scales within galaxies, especially at low 
surface densities. 

Estimates of the slope of the SK relation within individual 
galaxies exhib i ts sign ifica nt variations. For exa mple, while 
ISchuster et all (120071) and iKennicutt etall (120071) find n tot m 
1.4 for the molecular-rich galaxy M51a, similar estimates in 
other large, nearby ga laxies (including the Milky Way (MW); 
Misiri otis et alj 12006)) range from n tot m 1.2 to « tot « 3.5 
dWong & BiitJ2 002: Boi ssier et al.ll2003l) depending on dust- 
corrections and fitting methods. The disk-averaged total gas 
SK relation for normal (non-starburst) galaxies also has a 
co mparably steep s lope of n tot « 2.4, with significant scat- 
ter (lKennicutj [l998). While the variations in the EsFR-£ gas 
correlation may indicate systematic uncertainties in observa- 
tional measurements, intrinsic variations or trends in galaxy 
properties may also in duce differences b etween the global re- 
lation determined by iKennicutl ([1998) and the Ssfr _ E gas 
correlation in individual galaxies. 

Galaxies with low gas surface densities, like dwarfs or bul- 
geless spirals, displa y an even w i der va riatio n in their star 
forma tion relations. iHever et ail d2004l) and iBoissier et alj 
(2003) show that in low-mass galaxies the SFR dependence 
on the total gas surface density exhibits a power-law slope 
»fnt ~ 2 — 3 tha t is considerably steeper than the global 
IKennicutl d!998l) relation slope of n tot « 1 .4. Further, star 
formation in the low-surface density outskirts of galaxies also 
may not be universal. Average SFRs appear to dr op rapidly at 
gas surface densities of Sg as < 5 — 1 M pc~ 2 dHunter et alj 
119981: iMartin & KennicutfcOOll) . indicating that star forma- 
tion may be truncated or exhibit a steep dependence on the 
gas surface density. The existence of such threshold sur- 
face densities have been prop osed on theoretical grounds 
dKennicuttll 19981: ISchaydl2004i) . although recent GALEX re- 
sults using a UV indicator of star formation suggest that 
star formation may continue at even lower surface densities 
(Boissier et al. 2007). Star formation rates probed by damped 
Lyman alpha absorption ( DLA) systems also appear to lie be- 
low the [Kennicuttl d 19981) relation, by a n order of magnitude 
dWolfe & Chenll2006t IWild et alj|2007l) . which may indicate 
that the relation between SFR and gas surface density in DLA 
systems differs from the local relation measured at high S gas . 

In contrast, observations generally show that star forma- 
tion in galaxies correlates strongly with molecular gas, es- 
pecial ly with the highly dense gas traced b y HCN emis- 
sion dGao & Solomon! 12004 IWu et al l 120051) . The power- 



law index of the SK relation connecting the SFR to the sur- 
face density of molecular hydrogen consistently displays a 
value of n mo \ « 1.4 and exhibits considerably less galaxy- 
to-galaxy variation ( [Wong & BlitzH2002t iMurgia et alJ | 2002 ; 
Boissier et al. 2003; Heyer et al.l|20 04; Matt hews etaTlboOS : 
Leroy et al.ll2005L l2006fc IGardan et al]l2007l) . Molecular gas, 
in tur n, is expected to form in the high-pressure regions o f the 
ISM dElmegreenlll993al: lElmegreen & Parravanol ll994l), as 
indicated by obser vations dBlitz & Rosolowskv 2004, 120061 : 
IGardan et alJl2007l) . 

Analytical models and numerical simulations that tie 
star formation to the fraction of gas in the dense ISM 
are s uccessful in rep r oducing many observational trends 
(e.g. lElmegreenll2002l: [Kravtsovl 120031: iKrumholz & McKed 
[200l iLi et alj l2005bf p06t iTasker & Brvanl I2006L 120071: 
Krumholz & Thompson 2007; IWada & Normanll2007t iTassisI 
2007). Recently, several studies have explored star formation 
reci pes based on molec ular hydrogen. Pelup essv et alj d2006l) 
and iBooth et alj d2007l) implemented models for H2 forma- 
tion in gaseous disks and used them to study the molecular 
content and star formation in galaxies. However, these stud- 
ies focused on the evolution of galaxies of a single mass and 
did not address the origin of the SK relation, its dependence 
on galaxy mass or structure, or its connection to trends in the 
local molecular fraction. 

Our study examines the SK relation critically, including its 
dependence on the structural and ISM properties of galaxies 
of different masses, to explain the observed deviations from 
the global SK relation, to investigate other connections be- 
tween star formation and disk galaxy properties such as rota- 
tion or gravitational instability, and to explore how the tem- 
perature and density structure of the ISM pertains to the star 
formation attributes of galaxies. 

To these ends, we develop a model for the ISM and star for- 
mation whose key premise is that star formation on the scales 
of molecular clouds (^10 pc) is a function of molecular hy- 
drogen density with a u niversal SFE per free-fall time (e.g., 
Krumholz & Tan 2007). Molecular hydrogen, which we as- 
sume to be a proxy for dense, star forming gas, is accounted 
for by calculating the local H2 fraction of gas as a function 
of density, temperatu re, and metallicity using the photoion- 
ization code Cloudy dFerland et al.ll 1 998b to incorporate H2- 
destruction by the UV radiation of local young stellar pop- 
ulations. We devise a numerical implementation of the star 
formation and ISM model, and perform hydrodynamical sim- 
ulations to study the role of molecular gas in shaping global 
star formation relations of self-consistent galaxy models over 
a representative mass range. 

The results of our study show that many of the observed 
global star formation correlations and trends can be under- 
stood in terms of the dependence of molecular hydrogen abun- 
dance on the local gas volume density. We show that the 
physics controlling the abundance of molecular hydrogen and 
its destruction by the interstellar radiation field (ISRF) play a 
key role in shaping these correlations, in agreement with ear- 
lier calculations based on more idealized models of the ISM 
(Elmegreen 1993a; Elmeg reen & Parravanoll9 94). While our 
simulations focus on the connection between the molecular 
ISM phase and star formation on galactic scales, the for- 
mation of molecular hydrogen has also been recently stud- 
ied in simulations of the ISM on smaller, su bgalactic scales 
dGlover & Mac Low! 120071: iDobbs & Bonnelll 120071) . These 
simulations are complementary to the calculations presented 
in our study and could be used as input to improve the molec- 
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ular ISM model we present. 

The paper is organized as follows. The simulation method- 
ology, including our numerical models for the ISM, interstel- 
lar radiation field, and simulated galaxies, is presented in §|2] 
The results of the simulations are presented in $3] where the 
simulated star formation relations in galactic disks and corre- 
lations of the molecular fraction with the structure of the ISM 
are examined. We discuss our results in §j4]and conclude with 
a summary in §]5] Details of our tests of numerical fragmenta- 
tion in disk simulations and calculations of the model scaling 
between star formation, gas density, and orbital frequency are 
presented in the Appendices. Throughout, we work in the 
context of a dark-energy dominated cold dark matter cosmol- 
ogy with a Hubble constant Ho ~ 70km s _I Mpc _1 . 

2. METHODOLOGY 

To examine molecular gas and star formation in disk galax- 
ies, we develop a model that accounts for low temperature 
coolants (T < 10 4 K), calculates the equilibrium abundance 
of H2 in dense gas, and estimates the local SFR based on the 
local H2 density with a universal SFE per free fall time. This 
ISM model is detailed below in §12.11 and is referred to as 
the H2D-SF model (for "H2 density-star formation"). We ex- 
tend the H2D-SF model to account for the destruction of H2 
by an interstellar radiation field powered by local star forma- 
tion, and refer to the extended model as the H2D-SF+ISRF 
model (for "H2 density-star formation plus interstellar radia- 
tion field"). The H2D-SF+ISRF model is also detailed below 
in §12.11 The new molecular gas ISM and star formation mod- 
els are contrasted against a simple ISM model that includes 
only atomic cooling with a temperature floor T = 10 4 K and 
a SFR calculated from the total gas density above a density 
threshold, which we refer to as the GD-SF model (for "gas 
density-star formation"). The GD-SF model has been com- 
monly used in previous galaxy formation simulations. Each 
star formation and ISM model is numerically implemented in 
the smoothe d particle hydrodynamics (SPH) / N- body code 
GADGET2 (ISpringel et alJl200lt Upringel 2005), which is 
used to perform the simulations presented in this work. Gen- 
eral issues of numerical fragmentation in the disk simu lations 
of the kind presented here are described in §12.31 §12.41 and in 
the Appendix. The ISM models are applied to simulations 
of the isolated evolution of disk galaxies with a representa- 
tive range of circular velocities and structural properties. The 
galaxy models are designed as analogues of the well-studied 
nearby galaxies DD0154, M33, and NGC 4501, and are de- 
tailed in §12.61 The results of the simulations are presented in 

m 

2.1. ISM Properties 

The thermal properties of the ISM are largely determined 
by the radiative heating H and cooling A rates of interstellar 
gas. Given Ti and A, the rate of change of the internal energy 
u of gas with a density p g evolves as 



with additional terms to describe the energy input from feed- 
back or hydrodynamical interactions. The temperature depen- 
dence of Ti and A will depend on the metallicity and density 
of gas, and on the nature of the external radiation field. 

To calculate the above heating and cooling processes, as 
well as additional cooling processes from molecules and 
ionic species other than H2, we use the photoionization 



code C loudy version C06.02a, last described bv lFerland et al.l 
(1998). Cloudy calculates heating and cooling processes for 
H, He, and metals, as well as absorption by dust grains, pho- 
toelectric absorption, and molecular photoabsorption. The 
chemical network and cooling includes molecules and atomic 
ions that allow gas to cool to temperatures T < 10 4 K. For 
details on the molecula r network and cooling , the reader i s 
encouraged to exa mine EolEnbach & McKeg (fl97l ll989). 
iFerland et alJ (119981) . and the Cloudy documentation suite 
Hazy 4 . For purposes of calculating the efficiency of molec- 
ular line cooling using Cloudy, only the thermal line width 
is considered. At temperatures below our minimum temper- 
ature of (T ~ 100K), the influence of nonthermal motions 
on molecular line cooling will become increasingly important 
and should be considered when modeling lower temperature 
gas. 

For the heating term Ti, our Cloudy calculations include 
heatin g from the cosmic ultrav iolet background at z = fol- 
lowing lHaardt & Madaul (119961) . extended to include the con- 
tribution to the UV background from galaxies (see the Cloudy 
documentation for details). Cosmic ray heating and ioni zation 
are incorporated following lFerland & Mushotzky ( 1984]), and 
are important at densities where the ISM becomes optically- 
thick to UV radiation. The presence of a blackbody cos- 
mic microwave background radiation with a temperature of 
T w 2.7K is included. Milky Way-like (R v = 3.1) dust, in- 
cluding both graphite and silicate dust, is modeled with the 
local ISM abundance. Grain heating and cooling mech- 
anisms are included, follow ing Ivan Hoof et al.l (1200 ll) and 
IWeingartner & Draind (1200 ll) . 

The H2D-SF model of the ISM includes all the cooling 
and heating processes described above. In addition, the H2D- 
SF+ISRF model includes the following treatment of the ISRF, 
and is used to examine the consequences of soft UV radiation 
fr om star formation f or the molecular phase of the ISM. 

Mat his et alJ ([1983) modeled the ISRF needed to power the 
emission spectrum of dust both in the diffuse ISM and in 
molecular clouds in the Milky Way. The short-wavelength 
(A < 0.3yum) ISRF spectral energy distribution (SED) was 
inferred to scale roughly exponentially with a scale length 
R„ ~ 4kpc in the MW disk (for a solar galactocentric radius of 
Rq = lOkpc), tracing the stellar population, while the strength 
of the long-wavelength ISRF SED varied non-trivially with 
radius owing to the relative importance of emissio n from dust. 
The sh ort-wavelength ISRF SED determined by Math is et alJ 
d 19831) (their Table A3) for R = 0.5- 1.3# , normalized to 
a fixed stellar surface density, is roughly independent of ra- 
dius, w i th a sh ape similar to the radia t ion fie ld considered by 
iDraind d!978l) and lDraine & Bertoldil d 19961) . For simplicity, 
we fix the ISRF spectrum to have its inferred SED in the so- 
lar vicinity and scale the intensity with the local SFR density 
normalized by the solar value as 

_ u u E SF r 

l/iarf = ^ (2) 

where SspR.rri » (2-5) x lO^M^yr'pc"^ (e.g. : ISmith et alJ 
[19781: iMiiler & Scald fl979t lTalbot|[T980t lRanalll987l Il99lF 
Kroupa Il995h . For our model, we adopt Ssfr,© = 4 x 
lO" 9 M yr _1 pc" 2 . When included, the ISRF is added to the 
input spectrum for the Cloudy calculations as an additional 
source of radiation. 

4 http://www.nublado.org 
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FIG. 1. — Cooling (A) and heating (K) rates for interstellar and intergalactic gas as a function of gas density («h), temperature (T), metall icity (Z), and 
interstellar radiation field (ISRF) strength (LW), in units of the ISRF strength in the MW at the solar circle, as calculated by the code Cloudy ( Ferl and et alj 
1998). Shown are the cooling (dotted lines), heating (dashed lines) and net cooling (solid) functions over the temperature range T = 10 2 -10 9 K. Dense gas can 
efficiently cool via atomic and molecular coolants below T = 10 4 K, depending on the gas density and the strength of the ISRF. A strong ISRF can enable the 
destruction of H2 gas and thereby reduce the SFR. 



Absorption of soft UV (A ~ 0. 1 /im) photons can enable the 
destruction of H2 through transitions to the vibrational con- 
tinuum or to excited states that can be p hotoionized or pho- 
todissociated (Stecher & Williams 1967). The ISRF can sup- 
ply the soft UV photons that lead to H2 dissociation and can 
subsequently regulate the H2 abundance at low gas densities 
(«h ~ 1cm" 3 ). For values of L W > 0.01, the ionizing flux 
of the ISRF dominates over the Haardt & Madau ( 1996) UV 
background (see, e.g. JSternberg et aLll2002f) . 

Figure Q] shows the cooling rate A (dotted line), the heat- 
ing rate TC (dashed line), and net cooling rate A ne t (solid 
lines) calculated for a range of temperatures (T), densi- 
ties («h), metallicities (Z), and ISRF strengths (f/i sr f). For 
each value of T, rtu, Z, and £/; sr f a Cloudy simulation is 
performed assuming a plane-parallel radiation field illumi- 
nating a lOpc thick slab chosen to be comparable to our 
hydrodynamical spatial resolution. Each simulation is al- 
lowed to iterate until convergence, after which the heat- 



ing rate, cooling rate, ionization fraction, molecular frac- 
tion, electron density, and molecular weight is recorded from 
the center of the slab. These quantities are tabulated on a 
grid over the range log 10 r = {2,9}log 1() K with Alog 10 r = 
0.25, log 1Q « H = {-6,6}log 10 cm -3 with Alog 10 n H = 1-5, 
\og m Z = {-2. l}log lo Z with Alog 10 Z = 1.5, and £/ isrf = 
{0, 1, 10, 100, 1000} t/i sl t,0, and log-linearly interpolated ac- 
cording to the local gas properties. We note that since the slab 
thickness determines the typical particle column density for a 
given volume density, the tabulated quantities can depend on 
the slab thickness. For reference, varying the slab thickness 
between lOOpc and lpc for a density of «h ~ 30cm" 3 changes 
the molecular fraction by less than 30% and the heating and 
cooling rates by less than 50%. 

The left panels of Figure [T] with Ui SI f = correspond to net 
cooling functions used in the H2D-SF model, while the H2D- 
SF+ISRF model includes all the net cooling functions pre- 
sented in the figure. The GD-SF model that includes only 
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atomic cooling incorporates the regions of the net cooling 
functions in the left (t/i sr f = 0) panels of Figure Q] at temper- 
atures T > 10 4 K. In addition to well known cooling and heat- 
ing processes operating at high temperatures (T > 10 4 K), the 
results of the Cloudy simulations show that cooling of gas 
near «h ~ lcm~ 3 is regulated by the presence of the ISRF 
(FigureQ] middle left and center panels) and an ISRF strength 
of t/i sr f ~ 1 contributes a low-temperature heating rate Ti more 
than an order of magnitude larger than that supplied by the 
cosmic UV background. At higher densities («h ~ 10 3 cm" 3 ), 
the gas remains optic ally- thin to cosmic ray heating but be- 
comes insensitive to the presence of either the ISRF or the 
cosmic UV background. The correlation between the ISRF 
field strength and the local SFR limits the applicability of the 
1ow-«h / high-f/i sr f region of the cooling function but we in- 
clude such regimes for completeness. 

The simple description of the ISRF used in our model- 
ing is designed to replicate average conditions in systems 
like the Milky Way. However, the ISRF in special loca- 
tions in the ISM, such as near photoionized regions, may 
have a different character than the average SED we employ. 
Also, changes in the composition of dust or the initial mass 
function of stars relative to Milky Way properties could al- 
ter the frequency-dependent ISRF. We therefore caution that 
the detailed, frequency-dependent connection between SFR 
and ISRF strength used in our modeling does not capture ev- 
ery condition within the ISM of a single galaxy or between 
galaxies. The model we present does realistically capture the 
physics that regulate the abundance of molecular hydrogen 
and local SFR for a given spectral form of the ISRF. 

2.2. Star Formation and Feedback 

Modeling the SFR in the ISM as a power-law function of 
the gas density p g dates to at least Schmidt (1959), who mod- 
eled the past SFR of the Galaxy needed to produce the ob- 
served luminosity function of main sequence stars. We as- 
sume that star formation occurs in molecular clouds in pro- 
portion to the local molecular gas d ensity fmPe, as suggested 



models is given by 
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by a variety of observat i ons (e.g. , Elmegreen & Lad; 
Blitz &Thaddeusl 1 19801; iBeichman et al l Il986t iLadal 
Young & Scovilldll991l) . The model thus assumes that H2 



1987 



is a good proxy for star forming gas. 

On the scale of individual gas particles the SFR is then de- 
termined by the density of molecular hydrogen, which is con- 
verted into stars on a time scale /*: oc fmPg/h- Observa- 
tions indicate that at high densities scales with the local free 
fall time of the gas, % oc p? 5 , as 



t&f £ff) 



(3) 



with the SFE per free fall time, eg » 0.02 , approximately inde- 
pendent of density at n 10 2 - 10 4 cm" 3 (iRrumholz & McKea 
l200HlKrumholz & Tanll2007h . 

We adopt U = 1 Gyr for gas at density «h = 10/i 2 cm" 3 
(% = 2.33 x 10 7 yrs), which corresponds to eg = 0.023. By 
calibrating the SFE to observations of dense molecular clouds 
through eff, our model differs from the usual approach of 
choosing gas consumption time scale to fit the global SFR and 
the SK relation in galaxies. As we show below, the chosen ef- 
ficiency results in global SFR for the entire galaxies consistent 
with observations (see §0. 

Incorporating the assumption that a mass fraction f3 of 
young stars promptly explode as supernovae, the SFR in our 



p\ = a-P)fmy-(-. 



«H 



0.5 



10/i 2 cm- 3 y ■ (4) 

We assume (3 ~ 0.1, appropriate for the Salpeter ( 1955) initial 
mass function. Implicit in Equation [4] is that the molecular 
fraction / H2 ma y vary with a variety of local ISM properties 

as 

fm = fm (p g , T, Z g ,U isT i) ■ (5) 

The H2D-SF model includes the dependence on the gas den- 
sity, temperature T, metallicity Z g . The H2D-SF+ISRF model 
additionally includes the dependence of the molecular frac- 
tion on the strength of the interstellar radiation field t/j sr f pa- 
rameterized as a fraction of the local interstellar field en- 
ergy density. In both the H 2 D-SF and H 2 D-SF+ISRF mod- 
els, the equilibrium molecular fractio n is t abulated by using 
the Cloudy calculations discussed in 32.11 Variations in the 
cosmic ray ionization rate or the dust model could affect the 
form of Equation [5] but we do not consider them here. For 
the simple GD-SF model that does not track the molecular 
gas abundance, we set fm = 1 in gas with densities above 
the star formation thr eshold density rin > 0.1cm" 3 (see, e.g., 
iGovernato et~alll2007l) . 

Figure [T] shows that the feature in the cooling rate at T ~ 
10 3 K for low metallicity, intermediate density gas (e.g., mid- 
dle left panel), which owes to molecular hydrogen, disappears 
in the presence of the interstellar radiation field. Equations 
Q] and |4] then imply that destruction of H 2 by soft UV pho- 
tons from young stars regulates star formation in the H 2 D- 
SF+ISRF model. The dissociation of H 2 by the ISRF means 
that ISM gas may be cold and dense but its SFR may be sup- 
pressed if the local ISRF is strong. We discuss this effect in 
the context of the simulated galaxy models in ^3] 

The energy deposition from supernovae into interstellar gas 
is treated as a thermal feedback, given by 



du 



(6) 



where £sn is the energy per unit mass of formed stars that 
is deposited into the nearby ISM. We choose e = 1.4 x 
10 49 erg s Mq or, in the language of Springel & Hernquisi 
(2003b), an effective supernovae temperature of 7sn = 3 x 
10 s . This value for the energy deposition for supernovae has 
been used repeatedly i n simulations of galaxy mergers (e.g., 
iRobertson et al.l 2006a b e). Given that extremely short cool- 
ing time of dense gas at low resolutions, prescribing super- 
novae feedback as thermal input into the gas has long been 
known to have a wea k effect on the global evolution of sim - 
ulated galaxies (e.g., lKatzHl992t ISteinmetz & Mulled [l995h . 
For this reason, the results of this paper also do not depend 
strongly on the value of £sn or the supernovae mass fraction 
(3. Metal enrichment of the ISM is treated in the instanta- 
neous approximation, with a mass fraction y = 0.02 of the 
supernovae ejecta being retu rned into the gas as metals (see 
ISpringel & Hernquistll2003al) . 

2.3. Avoiding Numerical Jeans Fragmentation 

Perturbations in a self-gravitating medium of mean density 
p and sound speed c s can g row only if their wavelength ex- 
ceeds the Jeans length Aj ean s dJeansf l928). The corresponding 
Jeans mass contained within Aj ea ns is 



^Jeans — 



6G 3/ V/2 



(7) 
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(see, e.g., §5.1.3 of iBinnev & Tremainel 1 1987b . 
iBate & Burkertl (1 19971) identified a resolution requirement for 
smoothed particle hydrodynamics simulations such that the 
number of SPH neighbors A^ ne i g h and gas particle mass m gas 
should satisfy 2/V ne i g h'M g as < '«jeans to capture the pressure 
forces on the Jean s scale and avoid numerical fragmentation. 
iKlein et ail (|2004|) suggest that the effective resolution of 
SPH simulations is the number of SPH smoothing lengths per 
Jeans scale 

7r 5 / 2 c 3 



6G 3 / 2 A^ n eigh«i gas/ 9 



1/2 > 



(8) 



which is very similar to the number of Jeans masses per SPH 
kernel mass, if we define the kernel mass as 



WSPH 



(9) 



Our simulations include low temperature coolants that allow 
ISM gas to become cold and dense, and for the typical number 
of gas particles used in our models (A^ gas s» 400, 000) the Jeans 
mass is not resolved at the lowest temperatures and largest 
galaxy masses. Motivated by techniques used to avoid numer- 
ical J eans fragmentation in grid codes (e.g., Ma chacek et al.l 
1200 lb . a density-dependent pressure floor is introduced into 
our SPH calculations. For every gas particle, the kernel mass 
is monitored to ensure it resolves some number Nj edns of Jeans 
masses. If not, the particle internal energy is altered following 



2/3 



Cleans < N] 



Jeans 



U -U X \ - 2/3 



(10) 



TWjeans < 2A/j ea ns'«SPH 



to provide the largest local pressure and to assure the local 
Jeans mass is resolved. We discuss the purpose and effect of 
this pressurization in more detail in the first Appendix, but we 
note that this requirement scales with the resolution and natu- 
rally allows for better-resolved SPH simulations to follow in- 
creasingly lower temperature gas. For the simulations results 
presented in this paper, we use iVjea ns =15, which i s large r 
than the effective /V Je ans = 1 used bv lBate & Burkertl d 19971) . 
We have experimented with simulations with Nj eans = 1-100 
and find Afj ean s ~ 1 5 to provide sufficient stability over the time 
evolution of the simulations for the structure of galaxy mod- 
els we use. Other simulations may require different values 

-2/3 

Njeaas or a different pressurization than the u oc wi Jeans scaling 
suggested by Equation|7] 

2.4. Median Equation-of-State 

In the absence of a numerical pressure floor, the equation- 
of-state (EOS) of gas will follow the loci where heating bal- 
ances cooling in the temperature-density phase space. How- 
ever, the gas pressurization requirements to avoid numerical 
Jeans fragmentation change the EOS of the gas. For purposes 
of calculating the molecular or ionization fraction of dense 
ISM gas that has been numerically pressurized, an estimate 
of the EOS in the absence of resolution restraints is required. 
Figure|2]shows the temperature-density phase diagram of ISM 
gas in a simulation of a Milky Way-sized galaxy without us- 
ing the SPH pressure floor described in ^2.31 and detailed in 
the Appendix. 

The median EOS of this gas, determined by the heating 
and cooling processes calculated by the Cloudy code for a 
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FIG. 2. — Median equation-of-state (EOS) of gas in a Milky Way-sized 
galactic disk in our ISM model. Shown is the histogram of particles in the 
hydrogen number density nu - temperature T plane and the median temper- 
ature measured in bins of density (white line). The median EOS is used to 
determine the p rope rties of g as pressurized to avoid numerical Jeans frag- 
mentation (see j|2.3l and j|2.4t . At high densities, the particles follow the loci 
where heating balances cooling for their densities, temperatures, metallici- 
ties, and local interstellar radiation field strengths. 

given density, temperature, metallicity, and interstellar radia- 
tion field strength, as well as heating from supernovae feed- 
back, is used in our simulations to assign an effective tem- 
perature to dense gas pressurized at the Jeans resolution limit. 
The behavior of gas above a density of «h ~ 100cm" 3 is es- 
pecially influenced by the balance of available cooling mech- 
anisms and supernovae heating, and in the absence of super- 
novae feedback such dense gas would cool to the minimum 
temperature treated in our simulations (T ~ 100K). 

2.5. Molecular Fraction vs. Gas Density 

The molecular content of the ISM in models that account 
for the abundance of H2 and other molecular gas depends on 
the physical properties of interstellar gas through Equation 
[5] Since the temperature and density phase structure of the 
ISM is dictated by the local cooling time, which is nonlin- 
early dependent on the temperature, density, metallicity, and 
ISRF strength, the calculation of the equilibrium molecular 
fraction as a function of the gas properties is efficiently cal- 
culated using hydrodynamical simulations of galaxy evolu- 
tion. We have measured the molecular gas fraction fm as 
a function of gas density «h after evolving in ISM models 
with H 2 -destruction by an ISRF (the H 2 D-SF+ISRF model), 
and without an ISRF (the H 2 D-SF model), for t = 0.3Gyr 
for each of the v c j rc = 50- 300 km s" 1 disk galaxy simula- 
tions. The galaxy-to-galaxy and point-to-point variations in 
the metallicity, temperature, ISRF strength contribute to the 
spread in molecular fractions at a given density «h- With- 
out H2-destruction, the molecular hydrogen content of low- 
density gas calculated using Cloudy is surprisingly large, with 
fm = 0.1 at «h = 0.1cm" 3 and fm ~ 0.35 at nu = 1cm" 3 . For 
the H2D-SF model, the atomic-to-molecular transition repre- 
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sents a physical calculation of the star formation threshold 
density prescription often included in coarse models of star 
formation in simulations of galaxy formation and evolution. 

For the H 2 D-SF+ISRF model, the destruction of H 2 by 
the ISRF reduces the molecular fraction in the low-density 
(« H < 1cm" 3 ) ISM to / H 2 < 0.05. The abundance of H 2 in 
this model is therefore greatly reduced relative to the H2D-SF 
model that does not include an ISRF. The decline in molec- 
ular abundance, in turn, produces much lower SFRs in low 
density gas (rin < 1cm" 3 ) compared with either the H2D-SF 
model (which has a larger molecular fraction at low densities) 
or the GD-SF model (which allows for star formation in all 
gas above a density threshold, «h > 0.1cm" 3 ). 

For denser gas in the H2D-SF+ISRF model (e.g, «h ^ 
lcm" 3 ), the numerical pressurization of the ISM required to 
mitigate artificial Jeans fragmentation prevents much of the 
gas from reaching high densities where the molecular frac- 
tion becomes large (fm ~ 1 at «h ^ 100cm" 3 ). The molec- 
ular fraction of the numerically-pressurized gas must be as- 
signed at densities that exceed resolution limit of the simula- 
tions. As Figure [^indicates, the transition between the warm 
phase at densities «h < lcm" 3 a nd the cold phase at densities 
«h ^ 30cm~ 3 is rapid (see also IWolfire e t al. 2003). The re- 
quired pressure floor can artificially limit the density of gas to 
the range «h ^ 30cm" 3 where, in the presence of an ISRF, the 
molecular fraction of the gas would be artificially suppressed 
(without an ISRF, even gas at «h ~ lcm" 3 is already mostly 
molecular and would be relatively unaffected). To account for 
this numerical limitation, gas at densities «h > lcm" 3 that is 
numerically-pressurized to resolve the Jeans scale is assigned 
a minimum molecular fraction according to the median EOS 
and the high-density («h > 30cm" 3 ) /h2 - «h trend of gas in 
the absence of the pressure floor. The fm - "h relation cal- 
culated by Cloudy can be modeled for pressurized gas in this 
regime as f m = max{[0.67 (log, n H + 1-5) - 1], 1}. 

The utilization of the median EOS to assign detailed prop- 
erties to the dense gas in this manner has the further benefit of 
adaptively scaling with the local resolution, allowing an ap- 
proximation of the ISM properties that improves with an in- 
crease in the particle number used in a simulation. Through- 
out the paper, when simulation results include gas particles 
that are pressurized to avoid numerical Jeans fragmentation, 
their reported temperatures and molecular fractions are deter- 
mined from the median EOS presented in Figure [2] 

While the density dependence of molecular fraction utilized 
here may not be correct in every detail (i.e., the transition to 
/h2 = 1 may occur at an inaccurately low density), for our pur- 
poses the primary requirement is that the model captures the 
fraction of star forming gas as a function of density in a realis- 
tic fashion. The gas mass in disks in the H2D-SF+ISRF model 
at E gas > 10M Q pc" 2 is mostly molecular and the molecu- 
lar fraction increases realistically with the external ISM pres- 
sure (see 83.31>. in good agreement with observations (e.g., 
IWong & BlitzEoollBlitz & Rosolowskvll2006l) . 

2.6. Isolated Galaxy Models 

We use three cosmologically-motivated galaxy models to 
study global star formation relations. Below we describe the 
models and list the observational data used to motivate their 
adopted structure. 

The gala xy initial conditions follow the methodology de- 
scribed by Springel et al. (2005), with a few exceptions. The 
properties of the disks are chosen to be representative of 
disk galaxies that form within the CDM structure formation 



paradig m (e.g.. | Mo et alj fl998l) . The model galaxies con- 
sist of a Hem quisddl990h dark matter halo, exponential stel- 
lar and gas eous disks with s cale lengths /? d ,* and R^, and 
an optional Hernquist (1990) ste llar bulge with s cale radius 
flb,+- For a given virial mass, the iHernquisll (Il990|) dark mat- 
ter hal o parameters are scaled to an effective iNavarro et alJ 
(1996) halo concentration. The velocity distributions of the 
dark matter halos are init ialized from the isotropic distribu- 
tion function provided in Hernquist (1990), using the rejec- 
tion technique described by iPress et al.l (119921) . The vertical 
distribution of the stellar disk is modeled with a sech 2 (z/2z d ) 
function with z d = 0.2 — 0.3i?d,*. The stellar disk velocity field 
is initialized by using a Gaussian with the local velocity dis- 
persion a\ determined from the potential and density via the 
Jeans equations. The stellar disk rotational velocity disper- 
sion and streaming velocity are set wit h the epicyclical ap- 
proxi mation (see §3.2.3 and §4.2. 1(c) of iBinney & Tremaina 
119871) . The bulge velocity field is also modeled by a Gaussian 
with a velocity dispersion determined from the Jeans equa- 
tions but with no net rotation. The gas disk is initialized as 
an isothermal medium of temperature of T = 10 4 K in hydro- 
static equilibrium with the total potential, including the disk 
self-gravity, which d etermines vertical d istribution of the gas 
self-consistently (see Spr ingel et alJl2005l) . The gravitational 
softening lengths for the dark matter and baryons are set to 
eoM = 100/i"'pc and e\, aryon = 50/r'pc, respectively. 

The galaxy models are designed to roughly match the ob- 
served rotation curves, and HI, H2, and stellar surface mass 
distributions of DD0154 (v circ w 50km s" 1 ), M33 (v circ w 
125 km s" 1 ), and NGC4501 (v drc « 300km s" 1 ). These sys- 
tems cover a wide range in circular velocity, gas fractions, 
disk scale lengths, molecular gas fractions, gas surface densi- 
ties, and gas volume densities. The parameters of the galaxy 
models are provided in Table Q] 

2.6.1. DDO 154 Analogue 

DDO 154 is one of the most gas rich systems known and 
therefore provides unique challenges to models of the ISM 
in galaxies. For the DDO 15 4 galaxy analogue, we used the 
Carignan & Beaulieul (119891) HI map and rotation curves as 
the primary constraint on the ma ss distribution. The total gas 
mass is M d . g sa 2.7 x lQfy^ Afo dCarignan & Beaulieu|[l989l: 
iHunter & Elmegreenll2004h . The stel lar disk mass h as been 
estimated at M dj * w3.4x lO 6 h~ l M dLee et alJ l2006) with a 
scale length of about R±+ w 0.38/i~'kpc at A = 3.6/im. The 
gaseous disk scale length was determined by approximating 
the HI surface density profile, with Rd, g ~ 1.52"'kpc provid- 
ing a decent mimic of t he HI data. These numbers are simi- 
lar to those compiled bv lMcGaughl d2005l) for DDO 154 and, 
combined with the dark matter virial velocity V200 = 50 km s" 1 
and concentration c = 6, produce a rotational velocity of 
Vdrc ~ 50km s" 1 (where v c ; lc is defined as the maximum of ro- 
tation veloc ity profile) similar to the va lue of v c j rc as 54 km s" 1 
reported by lKarachentsev etlill (l2004h . 

2.6.2. M33 Analogue 

The nearby spiral M33 has well measured stellar, HI, H2, 
and SFR distributions, and serves nicely as an example disk 
galaxy with an intermediate r otational velocity. F or the M33 
galaxy analogue, we used the iHeyer et all (12004)) HI and H2 
map as guidance for determining the gas dis tribution with a 
total gas mass of M dig « 1.68 x 1O 9 A _1 M dCorbellil 12001 
iHever etaT]|2004h . According to this data, the gas disk scale 
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TABLE 1 
Galaxy Models 



Structural Parameters 

"Galaxy P200 Vd^ c M disk 7g^ Wl Wg Mbulge «M 

Analogue km s km s~' /i _1 Mg /r'kpc /i~'kpc Ii^'Mq ft~'kpc 

DDO 154 50 50 6 2.91 X 10 s 0.99 0.38 1.52 

M33 100 125 10 4.28 X 10 9 0.40 0.98 2.94 

NGC4501 180 300 14 1.02x 10" 0.04 3.09 2.16 1.14xl0 10 0.62 

Numerical Parameters 
Ualaxy Af^J 77J7, N bulgt , 

Analogue /i~'kpc /r'kpc h Mq Data Refs. 

DDO 154 120000 20000 400000 ~ oTl O05 64 4.60 X 10" 1,2,3,4,5 

M33 120000 120000 400000 ■■■ 0.1 0.05 64 2.74 X 10 5 5,6,7 

NGC4501 120000 177600 400000 22400 0.1 0.05 64 6.51 x 10 5 8,9,10,11,12 



Ref erenc es. — ( 1) Carianan & Beaulieu (1989); (2) iHunter & Elmeereen (2004); (3) Hunter et al] 
(200 6ft: (4) iLee et al.1 j2006h; (5) iMcGaughl d200g); (6) iCorbellil (120031); (7) [Hever et alj d2004l): (8) 
IWong & Blitzl d2002l); (9) Mollenho ff & HeidtM200l|) ( 1 OnGuhathakurt a etal] ( 119881) ; (11) iRubin et alj 
dl999f) ; (12) lBoissier et all d2003l) 



length is roughly R^ g ~ 2.7/T 1 kpc. We set the total stellar 
disk mass to M a .* ~ 2.6 x 10 9 /z _1 M Q to match the total bary- 
onic mass of Mb ary on ~ 4.3 x l0 9 h~ l M Q reported by McGaugh 
d2005l) . and set the disk scale length to i? d ,* « 0.9/r'kpc. 
These numbers are similar to those compiled by McGaugh 
d2005l) for M33 and, combined with the dark matter virial ve- 
locity V200 = 100 km s -1 and concentration c = 10, produce a 
rotational velocity of v c i rc w 125 km s . 

2.6.3. NGC 4501 Analogue 

The massive spiral N GC 4501 is a very lumi nous disk 
galaxy (M K = -25.38, IMollenhoff & Heidtl 120011) that is 
gas p oor by roughly a factor of ~ 3 for its stellar 
mass dGiovanelli & Havnesl 119831; iKennev & Youngl 119891; 
Blitz & Rosolowskv 2006), and serves as an extreme con- 
trast to the dwarf galaxy DDO 154. To construct an 
analogue to NGC 4501, we use d the rotation curve data 
compiled by iBoissier et al.l d2003l) from the obse rvations of 
iGuhafhakurta et all (1 1988b and lRubin et al.l d 19991) as a guide 
for designing the potential. The stellar disk scale length Rd * as 
3.1/z _1 kpc and bulge-to-disk ratio Mbui ge /Mdi 5 k = 0.112 were 
taken from the K—bund observations of IMollenhoff & Heidtl 
d2001l) . and th e total stellar mass w as chosen to be con- 
sistent with the iBell & de Jongl d2001l) K-band Tully-Fisher 
relation. The gas distribution was modeled after the sur - 
face mass density profile presented by Wong & Blitz (2002). 
These structure properties, combined with a virial velocity 
V200 = 180 km s" 1 and dark matter concentration c = 14, pro- 
vide a circular velocity of v c i rc w 300 km s" 1 appropriate for 
approximating the observed rotation curve. 

3. RESULTS 

The isolated evolution of each of the three disk galaxy mod- 
els is simulated for t = l.OGyr using the three separate ISM 
and star formation models described in §|2] the model H2D- 
SF with H 2 -based star formation, the model H 2 D-SF+ISRF 
which in addition includes the gas heating by local interstellar 
radiation, and the GD-SF model with a T = 10 4 K temperature 
floor and constant threshold density for star formation. In all 
cases, the ISM of galaxies is assumed to have initial metallic- 
ity of 10~ 2 Z Q and is self-consistently enriched by supernovae 
during the course of simulation. 

3.1. Galactic Structure and Evolution 



The gas distributions in the simulated disk galaxies after 
t = 0.3Gyr of evolution in isolation are shown in Figure [3] 
The image intensity reflects the disk surface density while 
the color indicates the local mass-weighted temperature of the 
ISM. A variety of structural properties of the galaxies are ev- 
ident from the images, including ISM morphologies qualita- 
tively similar to those observed in real galaxies. The large- 
scale structure of each disk model changes slowly over gi- 
gayear time-scales, excepting only the H2D-SF model for the 
Vcirc = 50km s" 1 galaxy. In this case, the ISM fragments ow- 
ing to the large cold gas reservoir that becomes almost fully 
molecu lar in the absence of H2-destruction (see 33.2.11 and 
33.2.2l for further discussion). 

The simulated disks can develop spiral structure and bar- 
like instabilities during their evolution, especially when mod- 
eled with T < 10 4 K coolants. The detailed structures of the 
spiral patterns and bar instabilities are likely influenced by 
numerical fluctuations in the potential, but they illustrate the 
role of spiral structure in the ISM of real galaxies. As the ISM 
density is larger in these structures, the equilibrium gas tem- 
perature is lower compared with the average disk properties. 
These structures also contain higher molecular fractions (see 
33.3b and will be sites of enhanced star formation in our ISM 
model, similar to real galaxies. 

The small-mass systems evolved with the GD-SF model 
are more diffuse and have less pronounced structure than 
those evolved with low-temperature coolants, owing to the 
higher (T > 10 4 K) ISM temperature (and hence pressure) 
in high-density regions. The diffuse morphologies of low- 
mass galaxies in the GD-SF and H2D-SF+ISRF models owe 
to the abundance of warm atomic gas with thermal energy 
comparable to the binding energy of the gas. The tendency 
of low-mass galaxy disks to be more diffuse and extended 
than more massive systems owin g to inefficient cooling w as 
shown by iRravtsov et al.l d2004l) and iTassis et all (120081) to 
successfully reproduce the faint end of the luminosity func- 
tion o f galactic satellit e s and scaling relations of dwarf galax- 
ies. iKaufmann et al.l d2007l) used hydrodynamical simula- 
tions of the dissipative formation of disks to demonstrate that 
the temperature floor of T = 10 4 K results in the equilib- 
rium disk scale height-to-scale length ratio in systems with 
Vcirc ~ 40 km s" 1 to be roughly three times larger than for 
galaxies with v c ; rc ~ 80 km s . The disk morphologies in our 
GD-SF runs are qualitatively consistent with their results. 
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Fig. 3. — Gas distribution in galaxies with circular velocities of v C kc = 300km s" (left column), v c j rc = 125km s" (middle column), and 
Voire = 50 km s~ (right column) for differing models of the ISM and star formation. Shown are an atomic cooling model with a T = 10 4 K 
temperature floor and star formation based on the gas density above a threshold of ;jh = 0.1cm -3 ("GD-SF", bottom row), an atomic+molecular 
cooling model with a T = 10 2 K temperature floor and star formation based on the molecular gas density ("FbD-SF", middle row), and an 
atomic+molecular cooling model with a T = 10 2 K temperature floor, star formation based on the molecular gas density, and H2 -destruction by 
an interstellar radiation field ("H2D-SF+ISRF", top row). The intensity of the images reflects the gas surface density and is color-coded by the 
local effective gas temperature. The diffuse morphologies of low-mass galaxies in the GD-SF and H2D-SF+ISRF models owe to the abundance 
of warm atomic gas with temperatures comparable to the velocity scale of the local potential. 



3.1.1. Molecular Fraction vs. Disk Radius 

The variations in the molecular fraction fm as a func- 
tion of the gas volume density «h discussed in 32.51 map 
into /h2 gradients with surface density or disk radius on the 
global scale of a galaxy model. Figure [4] shows azimuthally- 
averaged molecular fraction fm as a function of radius nor- 
malized to the disk scale length in the simulated disk galaxies 
with v c ; rc = 50- 300km s" 1 for molecular ISM models H2D- 
SF+ISRF and H 2 D-SF, after t = 0.3 Gyr of evolution. The av- 
erage molecular fraction decreases with decreasing rotational 
velocity, from fm > 0.6 within a disk scale length for the 
Vcirc = 300km s" 1 system (left panel) to fm < 0.4 within a 
disk scale length for the v a kc = 50 km s" 1 system (right col- 
umn) in the H2D-SF+ISRF model. The radial dependence of 
the molecular fraction steepens in smaller mass galaxies rel- 
ative to the Voire = 300km s" 1 model since a larger fraction 
of gas has densities below the atomic-to-molecular transition. 
The impact of this trend on the global star formation rela- 
tion is discussed below in 33.2.11 We note that details of the 



model initial conditions can influence the radial distribution 
of the molecular gas by changing the typical gas density in 
the disk. In particular, the v irc = 125 km s" 1 model galaxy has 
a more centrally concentrated molecular gas distribution than 
does M33. 

Figure g] clearly shows that the ISRF in the H 2 D-SF+ISRF 
model steepens the dependence of molecular fraction fm 
on the radius R for each galaxy compared with the H2D-SF 
model. 

3.1.2. Star Formation Histories 

The star formation histories of the simulated galaxies pro- 
vide a summary of the efficiency of gas consumption for 
systems with different disk structures and models for ISM 
physics. Figure [5] shows the SFR over the time interval t = 
0-1 Gyr for the isolated galaxy models evolved with the GD- 
SF, H 2 D-SF, and H 2 D-SF+ISRF models. In most simulations, 
the SFR experiences a general decline as the gas is slowly 
converted into stars. Substantial increases in the SFR with 
time are driven either by the formation of bar-like structures 
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FIG. 4. — Molecular fraction as a function of radius in three simulated disk galaxies with ISM models with (light gray lines) and without (dark gray lines) an 
interstellar radiation field. The general trend is for molecular fraction to decrease and its dependence on radius to steepen with decreasing mass of the system. 
This trend is more pronounced in the model with an interstellar radiation field. 
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FIG. 5. — Star formation rates (SFRs) as a function of time in isolated disks. Shown are the SFRs for galaxy models with v c j rc = 50km s (right panel), 
Voire = 125 km s (middle panel), and v c j rc = 300km s (left panel) calculated from interstellar medium models that allow for cooling at temperatures T > 10 4 K 
(black lines) or cooling at temperatures T > 10 K with (light gray lines) and without (dark gray lines) destruction of H2 by a local interstellar radiation field 
(ISRF). The star formation efficiency of the low-mass systems is strongly influenced by H2-destruction from the ISRF. 

(e.g., after t as 0.7Gyr for the v c i rc = 300km s" 1 galaxy) or, in 
one case, fragmentation of the cold ISM (after t as 0.4Gyr for 
the Vcirc = 50km s" 1 galaxy evolved with the H2D-SF model). 

For the most massive galaxy, differences in the relative star 
formation efficiencies of models that tie the SFR to the total 
gas density (the GD-SF model) or to the molecular gas density 
(the H 2 D-SF and H 2 D-SF+ISRF models) are not large. The 
ISM in the massive galaxy is dense, and when H 2 abundance 
is modeled much of the gas becomes molecular (see Figured . 
When H2-destruction is included the molecular content of the 
ISM in this galaxy is reduced at large radii and small surface 
densities, but much of the mass of the ISM (e.g., interior to 
a disk scale length) retains a molecular fraction comparable 
to the ISM model without H 2 -destruction. The effect of the 
ISRF on the SFR of this galaxy is fairly small. 

The SFR history of the intermediate-mass (v c i rc = 
125 km s" 1 ) galaxy model is more strongly influenced by the 
H 2 -destroying ISRF. While the GD-SF and H 2 D-SF ISM 
models produce similar star formation efficiencies, reflecting 
the substantial molecular fraction in the ISM for the H 2 D- 
SF model, the effects of H 2 -destruction in the H 2 D-SF+ISRF 
model reduce the SFR by ~ 30-40%. The time derivative of 
the SFR in the H 2 D-SF+ISRF model is also shallower than for 



the GD-SF or H 2 D-SF models, reflecting a reservoir of non- 
star-forming gas. In each model, the SFR of the whole system 
is similar to that estimated for M3 3 (SFR~ 0.25-0.7 Af^yr" 1 
i Kennicuttl fl998t iHippelein et all l200l iBlitz & Rosolowskyl 
l2006tlGardan et al.ll2007l) . " 

The effects of H 2 -destruction play a crucial role in de- 
termining the star formation history of the v c ; lc = 50 km s -1 
dwarf system. The single phase ISM (GD-SF) model pro- 
duces a SFR history qualitatively similar to the more mas- 
sive galaxies, but in the H 2 D-SF model the star formation 
history changes dramatically. The colder ISM temperature 
allows the gaseous disk to become substantially thinner in 
this case compared to the GD-SF model, leading to higher 
densities and SFRs. After t as 0.4Gyr of evolution, the cold, 
dense gas in the H 2 D-SF model allows the ISM in the dwarf to 
undergo large-scale instabilities and, quickly thereafter, frag- 
mentation. The dense knots of gas rapidly increase their SFE 
compared with the more diffuse ISM of the GD-SF model 
and the SFR increases to a rate comparable to that calculated 
for the v circ = 125 km s" 1 model. In the H 2 D-SF+ISRF model, 
the equilibrium temperature of most diffuse (and metal poor) 
gas in the dwarf is increased and the vertical structure of the 
disk becomes thicker than in the H 2 D-SF model. The molec- 
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ular content of the disk is thus substantially reduced (see the 
right panel of Figure @]l and the disk remains diffuse. As a 
result, the SFR drops to a low and roughly constant level. The 
constancy of the SFR in the H2D-SF+ISRF model is a result 
of the large gas reservoir that remains atomic and is there- 
fore unavailable for star formation. As the molecular gas, 
which comprises a small fraction of the total ISM in the dwarf 
(fm ^0.1 by mass), is gradually consumed it is continually 
refueled by the neutral gas reservoir. 

The star formation history of the dwarf galaxy system in 
the H2D-SF+ISRF model demonstrates that physics other 
than energetic supernova feedback can regulate star forma- 
tion and produce long gas consumption time scales in dwarf 
galaxies. Star formation rates in real galaxies similar to 
the Vcire = 50km s _I mo del fall in the wide range SFR ~ 
lO^-O.lMoyr" 1 (e.g.. iHunter & Elmegreer] 120041) . While 
the H2D-SF+ISRF simulation of the dwarf galaxy produces 
a steady SFR — 0.025M Q yr _I , we note that the galaxy DDO 
154 that served as a model for the v c h-c = 50km s" 1 sys- 
tem has a star formation ra te of only SFR ~ 0.004M Q yr _1 
dHunter & Elmegreenll2004l) . 

3.2. Star Formation Relations in Galactic Disks 

The disk-avera ged SK relation in galaxies, determined by 
IKennicuttl d 1998b . is well described by the power-law 



£ SFR = (2.5±0.7)x 10" 



lM©pc _ 



1.4±0.15 



M Q yr kpc 



(11) 



However, as discussed in the Introduction, spatially- 
dependent determinations of the total gas SK relation slope 
vary in the range n tot = 1.7-3.55. While the slope of the 
spatially-resolved molecular gas SK relation is consistently 
measured to be n mo \ w 1.4 in the sa me systems, with a 
two-sigma variation of n mQ \ ~ 1.2 -1.7 dWong & Blitzl [20021 : 
Boissier et alj l2003h iHever et all [20041: iBoissier et alj 120071: 
Kenni cutt et alj2 007). Below, we examine SK relations in the 
simulated disks and compare our results directly with these 
observations. 

3.2. 1. The Total Gas Schmidt- Kennicutt Relation 

A dichotomy between the total gas and molecular gas SK 
relations is currently suggested by the data and should be a 
feature of models of the ISM and star formation. Figure [6] 
shows the simulated total gas SK relation that compares the 
SFR surface density Esfr to the total gas surface density £ gas 
for different galaxy models. The quantities Esfr and £ gas 
are azimuthally-averaged in annuli of width Ar = 150/i _1 pc 
for the v c „ c = 300 km s" 1 and v c j lc = 125 km s" 1 systems and 
Ar = 75/T 1 pc for the smaller v c irc = 50 km s" 1 disk. Regions 
within the radius enclosing 95% of the gas disk mass are 
shown. The SK relation is measured in each disk at times 
t = 0.3 Gyr and t = 1 .OGyr to monitor its long-term evolution. 
Power-law fits to the SFR density as a function of the total 
gas surface mass densities in the range £ eas = 5- 100M Q pc" 2 
measured at t = 0.3 Gyr are indicated in each panel, with the 
range chosen to approximate regimes that are currently obser- 
vationally accessible. The simulated relations between Ssfr 
and S gas are not strict power-laws, and different choices for 
the range of surface densities or simulation time of the fit may 
slightly change the inferred power-law index, with lower sur- 
face densities typically leading to larger indices and, similarly, 



later simulation times probing lower surface densities and cor- 
respondingly steeper indices. Exceptions are noted where ap- 
propriate, but the conclusions of this paper are not strongly 
sensitive to the fitting method of the presented power-law fits. 
Disk-averaged SFR and surface mass densities are also mea- 
sured and_gk2tted forcomparison with the non-starburst galax- 
ies of lKennicuttl (119981) . 

Each row of panels in Figure [6] shows results for the three 
galaxy models evolved with different models for the ISM 
and star formation. For the traditional GD-SF model, the 
Esfr _ S gas correlation in massive, thin-disk systems tracks 
the input relation for this model, p\ oc py , at intermediate 
surface densities. At low surface densities, massive disks de- 
viate from the Ssfr oc Si a 5 s relation as the imposed density 
threshold is increasingly probed in the disk interior and disk 
flaring begins to steepen the relation between the volume den- 
sity p g and the gas surface density £ gas . In lower mass galax- 
ies, the disks are thicker at a given ga s surface density, as dic - 
tated by hydrostatic equilibrium (e.g.. lKaufmann et al.ll2007l) . 
with a correspondingly lower average three-dimensional den- 
sity. The SK relation therefore steepens somewhat between 
galaxies with v c i rc = 300km s" 1 and v c i rc = 125 km s" 1 in this 
star formation model, but does not fully capture the mass- 
dependence of the SK relation apparent in the observations. 

The middle row of Figure|6]shows the total gas SK relation 
measured in disks evolved with the H2D-SF model, which are 
generally similar to the more prescriptive GD-SF model in 
terms of the relation slope. The arbitrary density threshold for 
star formation included in the GD-SF model roughly mimics 
the physics of the atomic-to-molecular gas transition calcu- 
lated by Cloudy for the H2D-SF model, as molecular transi- 
tion begins at densities of «h ~ 0.1cm" 3 when the effects of 
H2-destruction by an ISRF are not included. 

Note that the H2D-SF model allows gas to cool to consider- 
ably lower temperatures than in the atomic cooling-only GD- 
SF model. This allows the ISM gas in the v c i rc = 50km s -1 
galaxy to become cold enough to settle into a thin, dense disk 
and reach high molecular fractions. Indeed, the increase in 
the molecular content of the disk contributes to the increase in 
the SFE in the dwarf galaxy at t = 1 Gyr during its evolution 
seen in the middle right panel of Figure|6] However, observed 
Vcirc ~ 50km s" 1 galaxies have very low molecular fractions 
and the high molecular content of the dwarf at t = 1 Gyr is a 
potential deficit of the H2D-SF model. 

The bottom row of Figure [6] shows the SK relation for the 
most complete "H 2 D-SF+ISRF" model. In this model, the 
SK relation power-law index increases systematically from 
n tot = 2.2 for the v c i rc = 300km s" 1 galaxy to n tot = 4.3 for the 
v c i rc = 50km s" 1 dwarf. The decreased abundance of H2 at 
low densities has steepened the SK relation at low gas surface 
densities compared with either the GD-SF or H2D-SF mod- 
els. The ISM in the dwarf galaxy is prevented from becoming 
fully molecular and maintains an extremely steep SK relation 
over its evolution from t = 0.3 Gyr to t = 1 Gyr. Consequently, 
the global SFR of the model dwarf galaxy in this case stays 
low and approximately constant. The global gas consumption 
time scale is also much longer than in the other two models. 

While the disk-averaged SK relation determination for 
the galaxies have the lowest SFE in the H2D-SF+ISRF 
model, the values are still consisten t with the observationa l 
determinations bv IKennicuttl (1 19981). IWong & Blitzl (120021) . 
IBoissier et all (120031) . and lBoissier et alj (120071) . ~ 

3.2.2. The Molecular Gas Schmidt-Kennicutt Relation 
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Fig. 6. — Schmidt-Kennicutt relation for total gas surface densities. Shown is the SFR surface density Esfr as a function of the total gas 
surface density E gas for galaxy models with circular velocities of v c j lc = 50km s (right column), v c j lc = 125 km s (middle column), and 
Vote = 300 km s~ (left column) for ISM models that include atomic coolants (upper row), atomic+molecular coolants with (lower row) and 
without (middle row) Hi -destruction from a local interstellar radiation field. The surface densities Esfr and E gas are measured in a nnuli from 
the sim ulated disks after evolving for t = 0.3 Gyr (solid circles) and t = 1 Gyr (open circles) in isolation. The disk-averaged data from K ennicutil 
( 1998) is shown for comparison (small dots). The best power-law fits to the individual galaxies have indices in the range a ~ 1.7-4.3 (solid 
lines). Deviations from the Esfr oc E ga 5 s relation are caused by the radial variation in the molecular gas fraction /h2, the scale height of star- 
forming gasJisH^ndjhe scale height of the total gas distribution /i gas . The total disk averaged surface densities (grey diamonds) are consistent 
with the Kennicutt ( 1998) normalization for normal, massive star-forming galaxies (dashed lines). 

and the observational correlations we consider are azimuthal 
surface density averages in annuli. The two are not equiv- 
alent, as we discuss below in the next subsection (see also 
ISchave & Dalla Vecchiall2008l) . 

The disk-averaged SFR and molecular gas surface densi- 
ties in the H2D-SF model are consi stent with observ ations 
of the molecular gas SK relation by iKennicutH d 1998b . The 
molecular gas SK relation power-law index of this model is 
similar for galaxies of different mass n mo \ = 1.4—1.5. Un- 



As discussed in $3.21 the two-sigma range in observational 
estimates of the molecular gas SK relation power-law index 
is n mo \ « 1.2—1.7, and is smaller than the variation in the 
slope of the £sFR-£c, as relation. Figure [7] shows the molec- 
ular gas SK relation measured in galaxies evolved with the 
two molecular ISM models considered in this paper. The SFR 
and molecular gas surface densities are measured from the 
simulations as described in $3.2.11 and the reported power- 
law indices were fit to molecular gas surface densities in the 
range S H 2 = 0.1 - 100M Q pc" 2 . Note that the correlation of 
the SFR surface density with the H2 surface density is not an 
entirely trivial consequence of our H2-based implementation 
of star formation, as we implement our star formation model 
on the scale of individual gas particles in three dimensions 



like the Ssfr _ S gas relation, the molecular gas SK relation is 
well-described by a power-law over a wide range in surface 
densities in the H2D-SF model. As discussed previously, the 
low-mass dwarf galaxy experiences an increase in the SFE 
at a fixed surface density during the evolutionary period be- 
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FIG. 7. — Schmidt-Kennicutt relation for molecular gas surface densities. Shown is the star formation rate surface density Ssfr as a function of the molecular 
gas surface density Sh2 for galaxy models with circular velocities of v c j rc = 50km s -1 (right panels), i> c j rc = 125km s (middle panels), and v c ; rc = 300km s 
(left panels) and ISM models that include atomic+molecular coolants with (lower panels) and without (upper panels) Hi-destruction from a local interstellar 
radiation field. The surface densities Ssfr and £h2 are measured in annuli (dark circles) from the simulated disks after evolving for t = 0.3Gyr (solid circles) and 
t = lGyr (open circles) in isolation. The disk-averaged data from Kennicutt 1 1998) is shown for comparison (small dots). The best power-law fit s to the ind i vidual 
galaxies have indices in the range a tzi 1.2— 1.5 (solid lines). The total disk averaged surface densities (grey diamonds) are consistent with the Kennicutt 1 1998) 
relation for normal, massive star-forming galaxies (dashed lines). 



tween t = 0.3 Gyr and t = l.OGyr owing to a nearly complete 
conversion of its ISM to molecular gas in the H2D-SF model. 

Panels in the bottom row of Figure [7] show that the addi- 
tion of a local interstellar radiation field in the H2D-SF+ISRF 
model controls the amplitude of the Ssfr ~~ Sh2 relation, but 
does not strongly affect its slope. The normalization and 
power-law index of the molecular gas SK relation in this 
model are quite stable with time, reflecting how the equi- 
librium temperature and molecular fractions of the ISM in 
the H2D-SF+ISRF are maintained between simulation times 
t = 0.3Gyr and t = 1 .OGyr. Only the most massive galaxy ex- 
periences noticeable steepening of the Ssfr _ Sh2 correlation 
at high surface densities. The di sk- averaged obse rvations of 
the molecular gas SK relation by iKennicuttl (119981) are well- 
matched by the simulations and remain approximately con- 
stant over their duration. 

While the scale height of the cold, dense gas in the disks is 
typically comparable to our gravitational softening, the slope 
of the molecular gas SK relation appears to be insensitive 
to modest changes in the simulation resolution. Additional 
resimulations of the calculations presented here with an in- 
creased resolution of A^.g = 1.2 x 10 6 gas particles show that 
the slope of the molecular gas SK relation changes by less 
than 5%. 

3.2.3. Structural Schmidt-Kennicutt Relation 

While the observational total and molecular gas SK rela- 
tions are well-reproduced by the simulated galaxy models, lo- 



cal deviations from a power-law Ssfr _ S gas or Ssfr _ Sh2 
correlation are abundant. Although these deviations are physi- 
cal, exhibited by real galaxies, and are well-understood within 
the context of the simulations, a question remains as to 
whether a relation that trivially connects the SFR density Ssfr 
and the properties of the disks exists. 

The primary relation between the local star formation rate 
and the gas properties in 3D in the simulations are encapsu- 
lated by EquationHJ which relates the three-dimensional SFR 
density to the molecular gas density. We can relate the three- 
dimensional properties in the ISM to surface densities through 
the characteristic scale heights of the disks. The SFR surface 
density is related to an average local SFR volume density (p+) 
through the scale height of star-forming gas /isfr as 



^SFR 



oc (p'+) h 



SFR 



(12) 



Here, the averaging of the three-dimensional density is under- 
stood to operate over the region of the disk where the two- 
dimensional surface density is measured. Similarly, the gas 
surface and volume densities are related by the proportional- 
ity 



(Pg) h s 



(13) 



If Equation|4]is averaged over the scale on which Ssfr is mea 
sured as 

,1.5 



(p\) oc (/ H2 ) (p g ) 



(14) 
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FIG. 8. — A predicted star formation correlation for disk galaxies. The star formation rate (SFR) density in the simulations is calculated as p\ = fmPg/U oc 
fni(pg)Pg ' 5 > where /m is the molecular fraction, p g is the gas density, and U is a star formation time-scale that varies with the local gas dynamical time. The 
local SFR surface density should scale as Esfr oc p*/?sfr. where Asfr is the scale height of star-forming gas. Similarly, S ga s oc p g h g!ls , where /i gas is the total gas 
scale height. The simulations predict that the SFR density in disks should correlate with the gas surface density as Ssfr oc fmhsFRh^hfEb^ (dashed lines). 



then the proportionalities can be combined to give 

£SFR0C</H2)^£i' 5 



h 



1.5 gas- 



(15) 



Equation[l5]can be interpreted as a SK-like correlation, which 
depends both on the molecular gas abundance and the struc- 
ture of the disk. 

Figure [8] shows this expected structural SK relation for the 
simulated disk galaxies. The star formation, and surface den- 
sities are azimuthally-averaged in annuli, as before, and lim- 
ited to regions of the disk where star formation is active. The 
characteristic gas scale height /! gas is determined as a mass- 
weighted average of the absolute value of vertical displace- 
ment for the annular gas distributions. The characteristic star- 
forming gas scale height /isfr is determined in a similar man- 
ner, but for gas that satisfies the condition for star formation: 
namely, f m > for the H 2 D-SF and H 2 D-SF+ISRF models 
or «h > 0.1cm" 3 for the GD-SF model. Measurements from 
the simulations are presented after evolving the systems for 
t = 0.3Gyr with the GD-SF, H 2 D-SF, and H 2 D-SF+ISRF ISM 
models. In each case, as expected, Ssfr in the disks linearly 
con-elates with quantity / H2 /isFR/!g as 5 £g a 5 s 5 . 

This relation shows explicitly that the behavior of the 
Esfr - S gas relation can be understood in terms of the depen- 
dence of fui and /i gas on £ gas . For example, the SK relation 
can have a fixed power law £sfr oc £ ga 5 s wh en fm and /i gas 
are ro ughly independent of S gas (e.g., M5 1 a. iKennicutt et al.l 
120071) . but then steepen when the mass fraction of star form- 
ing gas becomes a strong function of gas surface density in 
a manner characteristic of dwarf galaxies, large low surface 
brightness galaxies, or the outskirts of normal galaxies. Devi- 
ations from the Eg a 5 s relation can also occur when /i gas changes 
quickly with decreasing E gas in the flaring outer regions of the 
disks. Such flaring, for e xample, is observe d for the Milky 
Way at R > 10 kpc (e.g.. IWolfire et al.ll2.QQ3L and references 
therein). 

3.3. Local Pressure, Stellar Radiation, and Wi-Fraction 

Elmegreen (1993b) calculated the dependence of the ISM 
molecular fraction / H2 on the local pressure and radiation field 

5 Small deviations in the GD-SF model at low surface densities owe to 
the sharp density threshold fo r star formation adopted in this model (see 
Schaye & Dalla Vecchia 2008, for a related discussion). 



strength, including the effects of H 2 self-shieldin g and dust 
extinction. The calculations of Elmegreen (1993b) suggested 
that 

(16) 



fm oc P^f 



where P ext is the external p ressure and / i s the v olume emis- 
sivity of the radiation field. IWong & Blitzl ([2002) used obser- 
vations to determine the correlation between molecular frac- 
tion fui and the hydrostatic equilibrium mid-plane pressure 
in galaxies, given by 



-GSL 



(17) 



where is the stella r velocity dispersi on and a is the gas ve- 
locity dispersion (see Elmegre enll989t) . From observations of 
seven disk galaxies and the Milky Way, Wong & Blitzl (120021) 
determined that the observed molecular fraction scaled with 



pressure as 



E h2 /Shi ocP' t 



,0.8 



(18) 

Subsequently, IBlitz & Rosolowskvl (120041) and 
iBlitz & Ros olowskv d2006l) st udied these correlations in 
a larger sample of 28 galaxies. iBlitz & Rosolowskvl d2006l) 
found that if the mid-plane pressure was defined as 



— = 272 cm -3 K 



km s 



M pc- 2 

-0.5 



M Q pc- 2 



0.5 



(19) 



where /i* is the stellar scale height, and if the local mass 
density was dominated by stars, then the molecular fraction 
scaled with 



WShiocP^ 2 . (20) 

IWong & Blitzl (120021) state that if the ISRF volume emis- 
sivity j oc X* oc S e and the st e llar vel ocity dispersion was 
constant, then the Elmegreen (1993b) calculations would 
provide a relation fm oc P hi much steeper than observed. 
IBlitz & Rosolowskvl (120061) suggested that if the stellar ve- 
locity dispersion scales as c* oc 5 and j oc £*, then the 
predicted scaling is closer to the observed scaling. These as- 
sumptions provide / H2 oc P e ' x 2 . 



Molecular Hydrogen and Star Formation 



15 



10' 



10' 



10 2 
10" 3 

10 2 
to' 



W 10 



10' 



10" 



I -~P U 



( 



H D-SF 

2 

• v =300 km/s 
iv'=125 km/s 
■ v° =50 km/s 



H 1 1 h 



H J- 



f ~P U 

H2 



H D-SF+ISRF 

a 

• v =300km/s 
av c =125 km/s 
■ v° =50 km/s 



10 



10' 



10' 



10' 



10° 



10" 



10' 



P / k [cm" K] 



Fig. 9. — Molecular fractions as a function of external pressure 
in isolated disks. Shown are the molecular-to-atomic gas surface 
density ratios Em/Em measured in ra dial bins as a function of the 
local external pressure _P cxt defined by Blitz & Rosolowsky d2004l) 
for galaxies with v c ir C ~ 300kms~' (circles), v c irc ~ 125kms~' 
(triangles), and v c i rc ~ 50 km s~ (squares) circular velocities. The 
simulations that include IF-destruction by an ISRF (lower panel) 
follow the observed Em/Em Pen power-law scaling (dashed 
line), while the simulations that do not include an ISRF (upper 
panel) have a weaker scaling than that found observationally 
(E H2 /S H i ~ P° x f, dotted line). 



The relation between £ H 2/£hi and P ext can be measured for 
the simulated disks. Figure|9]shows £h2/£hi as a function of 
P ext , measured from the simulations according to Equation[T9l 
in annuli, for simulated disks with (bottom panel) and without 
(top panel) ^-destruction by an ISRF. The simulations with 
an ISRF scale in a similar way to the systems observed by 
iBlitz & Rosolowsky! (120061) . while the simulated disks with- 
out an ISRF have a weaker scaling. The simulated galaxies 
have c* oc 5 and roughly constant gas velocity dispersions. 

Given the scaling of the ISRF strength in the simulations, 
for which £/; sr f oc £sfr oc £"g with n tot ^2-4 according to 
the total gas SK relations measured in 33.2.11 and the scale- 
dependent SFE, which will instill the proportionality £* cx 



Eg as with (3 y 1 in general, the scaling between the molecu- 
lar fraction and pressure fm oc P®£ is not guaranteed to be 
satisfied in galaxies with greatly varying SK relat ions. How- 
ever, gi ven the relation fm oc P%£ j~ l calculated by Elmegreen 
(1993b), a general relation between the molecular fraction and 
interstellar pressure can be calculated in terms of the power- 
law indices « of the total gas SK relation and (3 that describes 
the relative distribution gas and stars in the disk. 

Assuming a steady state, the volume emissivity in the sim- 
ulated ISM then scales with the gas surface density as 



7 0cE/ isrf cx£"«. 



(21) 



Given the ISRF strength in the simulations and a proportion- 
ality of the stellar surface density with the total gas surface 
density of the form £* cx £g as , we find that, with the scaling 
between the molecular f raction, pressure, an d ISRF volume 
emissivity calculated by lElmegreenl (Il993bl) . the fm - £ 
relation should follow 



gas 



fm oc £ 



2.2(1+/3/2)-h 
gas 



(22) 



If the molecular fraction scales with the mid-plane pressure as 
fm oc P a , then the index a can be expressed in terms of the 
power-law dependences of the SK relation and the £* - £ 
correspondence as 



gas 



a = 



2.2(1 + P/2)-n 
1+/3/2 ' 



(23) 



The largest galaxy simulated, with a circular velocity of 
^circ ~ 300 km s" 1 , has a total gas SK relation with a scal- 
ing £sfr oc £g as (see Figure [6]) and a nearly linear scaling of 
the stellar and gas surface densities ((3 w 0.97). According to 
Equation|23] power-law index for this system is then expected 
to follow a = 0.85. At the opposite end of the mass spectrum, 
the v C i rc = 50km s" 1 dwarf galaxy, which obeys a very differ- 
ent SK relation (« to t ~ 4) and £* - £ gas scaling (J3 ~ 4), has 
an expected a i=s 0.87. Th ese a values are remarkably close 
to the relation observed by B litz & Rosolowsky! (l2006h . who 
find a w 0.92. That the ISRF model reproduces the observed 
correlation in galaxies of different masses then owes to the 
opposite compensating effects of the scaling of j with £ gas 
and the £ + -£ gas scaling. 

For the simulations that do not include an ISRF, the scal- 
ing between molecular fraction and pressure must be indepen- 
dent of the ISRF volume emissivity. We find that /h2 ~ £g a * 
in these galaxies, with significant scatter, which provides the 
much weaker relation £h2/£hi ~ P&a- The simulations with- 
out an ISRF have fm -Pext scalings similar to this estimate or 
intermediate between this estimate and the observed relation. 

The good agreement between the fm - Pext scaling in the 
simulations with an ISRF and the observations provide an 
a posteriori justification for the chosen scaling of the ISRF 
strength in the simulations (Equation |2), even as the ISRF 
strength was physically motivated by the generation of soft 
UV photons by young stellar populations. The robustness of 
the /H2 _ ^ext may be a consequence of the regulatory effects of 
H2-destruction by th e ISRF, which moti vated in part the orig- 
inal calculations by Elmegreen (1993b) , in concert with the 
influe nce of the external pressure (Elmeg reen & Parravanol 
1994). However, we note that systems with a steep total 
gas SK relation, such as the v c ; rc ~ 50 km s" 1 and v orrc ~ 
125 km s" 1 systems, do break from the fm~Pext relation in the 
very exterior of the disk where the molecular fraction quickly 
declines. Observations of M33, the real galaxy analogue to 
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the v c i rc ~ 125 km s simulated galaxy, indicate that the sys- 
tem may deviate from the observed fm - Pext relation deter- 
mine d at higher surface densities in the very exterior of the 
disk dGardan et alj2007h . The surface density of M33 in these 
exterior regions is domina ted by gas rather than stars, and the 
iBlitz & Rosolowskyl (|2006l) relation may not be expected to 
hold under such conditions. The high-su rface density regions 
of M3 3 do satisfy the fm~ Pext relation (Blit z & Rosolowskyl 
2006), as does the simulated M33 analogue. We delay a more 
thorough examination of the exterior disk regions until further 
work. 

Recently, Boo th et all d2007h presented a model of star for- 
mation based on a subgrid model for the formation of molec- 
ular clouds from thermally unstable gas in a multiphase ISM. 
In their model, molecular clouds are represented as ballistic 
particles which can coagulate through collisions. This is quite 
different from the treatment of ISM in our model, as we do not 
attempt to model molecular clouds on small scales via a sub- 
grid model, but calculate the equilibrium abundan ce of molec- 
ular h ydrogen using the local gas properties. Boo th et alJ 
d2007l) show that their simulatio ns also reproduce the molec- 
ular f raction-pressure relation of IBlitz & Rosolowskyl d2004l 
2006). However, the interpretation of this relation in their 
model must be different from our interpretation involving the 
relation between the gas surface density, the stellar surface 
density, and the ISRF strength, because their model does not 
include the dissociating effects of the interstellar radiation 
field. 

In addition, in contrast to our model, th e global SK re- 
lation in simulations of Booth et al. (2007) does not show 
a break or steepening down to very small surface densities 
(X w 10~ 2 M pc~ 2 ; see their Fig. 15). Given the differences, 
it will be interesting to compare results of different models of 
H2-based star formation in more detail in future studies. 

3.4. Star Formation and Angular Velocity 

The SFR surface density Ssfr has been shown to correlate 
with the product of the gas surface density £ gas and the disk 
angular frequency fi as 

XSFR — CsfiSgasfi (24) 

with the constant C S o = 0.017 (lKennicuttlll998l) and where 
the angular frequency is defined as 



R dR 



(25) 



(e.g., §3.2.3 of iBinnev & Tremaine 1987). Star forma- 
tion relations of this form have been forwarded to c onnect 
the SFR to the cloud-cloud collision time-scale (e.g.. IWysel 
fl986l IWvse & SiUd [19891 Hull l2000h . to tie th e gas con- 
sumption time-sc ale to the orbital time-scale (e.g.. lSilklll997l; 
Elme greenl 1 1 9971) . or to reflect th e correlation betw een the 
ISM density and the tidal density (Elmegreen 2002). These 
ideas posit that the Esfr _ £ gas f2 correlation provides a rea- 
sonable physical interpretation of the oc yfp scaling. 

Figure [10] shows the SFR surface density Esfr as a func- 
tion of the total gas surface density consumed per orbit 
Csn^gasfi and molecular gas surface density consumed per 
orbit C^qSh?.^, compa red to the observed Ssfr - £ gas f2 re- 
lation Kennicutt| d!998l) . The local values of £ gas and Sh2 
are measured in radial annuli from the gas properties, while 
ft is measured from the rotation curve at the same loca- 
tion. The SFR surface density scales as £ gas f2 at Esfr > 
5 x lO~ 3 M0yr _1 kpc~ 2 . The offset between the disk-averaged 



Esfr- S gas 51 relation found by Kennicuttl d 19981) and the sim- 
ulated relation may owe in part to the typically higher values 
of in the disk interior that can affect the normalization of 
the azimuthally-averaged Ssfr - S gas 57 relation. This offset is 
comparabl e to the observe d offset between the disk-averaged 
relation of lKennicuttl (119981) and t he spatially-reso l ved S sfr ~ 
S gas fi relation found in M51a by iKennicutt et all d2007l) . As 
in the previous sections, scaling Ssfr - £h2^ works through- 
out the disk, while the Ssfr _ £ gas f2 scaling steepens signifi- 
cantly at lower surface densities where the molecular fraction 
is declining. 

That the SFR in the simulated disks correlates well with 
the quantity £ gas f2, even though the SFR in the simulations 
is determined directly from the local molecular gas density 
and dynamical time without direct knowledge of il, may be 
surprising. The correlation suggests that the local gas den- 
sity should scale with 17, although the physical reason for 
such a correlation in our simulations is unclear. While pre- 
vious examinations of the Esfr - £ gas £! correlation have de- 
veloped explanations of why 17 oc <J()d, these simulations have 
a coarser treatment of the ISM than is invoked in such mod- 
els; the molecular fraction and SFR are primarily functions 
of the local density. The 17 ~ ^fp2 correlation should then 
reflect the structural properties galaxy models rather than the 
detailed properties of the ISM that are somehow influenced 
by global processes in the disk. 

The galaxy models used in the simulations remain roughly 
axisymmetric after their evolution, modulo spiral structure 
and disk instabilities. Hence, the fundamental equation that 
describes the gravitational potential $ of the galaxies is the 
Poisson equation 



d 2 $ 
dz 



„ , d 2 $ 



(26) 



where the density p tota ] reflects the total density of the multi- 
component system. If the form of the potential generated by 
ptotai has the property that S7 2 = (l/R)d$/dR oc pa, then the 
quantity £ gas 17 will mimic the scaling of Ssfr- In the second 
Appendix, we examine the solutions for 17 determined from 
Equation [26] for the limiting cases of locally disk-dominated 
(ptotai « Pdisk) and halo-dominated (p total ps p DM ) potentials, 
and show that typically 17 oc ^fpdB(pd) where B(pd) ~ 0(1) 
is a weak function of density if the disk density is exponen- 
tial with radius. This behavior suggests that whether the disk 
dominates the local potential, as is the case for massive galax- 
ies, or the dark matter halo dominates the potential as it does 
in low-mass dwarfs, the angular frequency scales as 57 oc ^fp~2 
throughout most of the disk. Based on these calculations, we 
conclude that the correlation Esfr oc £ gas 17 in the simulations 
is likely a consequence of the density-dependent star forma- 
tion relation (Equation |4]i and the galaxy structure inducing 
the correlation 17 oc Jp~d. This correlation may be established 
during the formation of the exponential disk in the extended 
dark matter halo, and the correlation observed in real galax- 
ies therefore may have a similar origin. Star formation may 
then correlate with £ gas S7, but not be directly controlled or 
influenced by the global kinematics of the disk. 

3.5. Star Formation and Critical Surface Density Thresholds 

Disk surface density, dynamical instabilities, and the for- 
mation of molecular clouds and stars may be connected. A 
single-component thin disk will experience growing axisym- 
metric instabilities if its surface mass density exceeds the crit- 
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FIG. 10. — Star formation rate surface density Esfr as a function of the total gas surface density consumed per orbit CsnS gas f2 (left panel) and molecular gas 
su rface density cons umed per orbit CsqEh2^ (right panel), where Q is the orbital frequency and the efficiency Cgn = 0.017 was determined from observations 
by Kennicutt 1 1998) (represented by the dashed line which extends over the range probed by observations). Shown are the local values of E gas , Sh2 and f! 
in annuli for galaxy models with circular velocities of v c j rc = 50km s (squares), v t .j rc = 125km s (triangles), and v c j rc = 300km s -1 (circles) evolved with 
atomic+molecular cooling and destruction of H2 by a local interstellar radiation field for t = 0.3Gyr. The star formation rate in the simulations display the 
correlations Esfr oc Sh2^ and, in the high surface density region where the molecular gas and total gas densities are comparable, Esfr oc Egasfi. 



ical value Eg, given by 

Efl(o-) : 



where k is the epicyclical frequency 

- 2 v dv 



2 -r+ 

R 2 RdR 



(27) 



(28) 



a is the characteristic velocity dispersion of the fluid, and the 
parameter a e ~ 0(1) (ISafronovlll9"60t lToomrelll964 . For a 
two component system, consisting of stars with surface mass 
density and velocity dispersion c*, and gas with surface 
mass density S gas and velocity dispersion a, the instability 
criterion for axisymmetric modes with a wavenumber k can be 
expressed in terms of the normalized wavenumber q = kc*/K 
and the ratio of velocity dispersions R = a/c+ as 



LI 



2 

Qs l+q i + Qz ~l+q 2 R 2 



> 1 



(29) 



dJog & Solomonll984tlRafikovl2001h . where Q s = S e (c*)/S i 
and Q g = Sg(cr)/Eg as . The most unstable wavenumber q can 
b e determined by m aximizing this relation. 

lElmegreenl(ll993al) considered another critical surface den- 
sity Sa, given by 

- - ° aA(7 (30) 



where the parameter 



irG 

0(1), and the Oort constant A is 



A = -0.5R- 



dR 



(31) 



(e.g.. lBinney & TremaindfT987h . above which density pertur- 
bations can grow through self-gravity before being disrupted 
by shear. This instability criterion ha s been shown to work 
well in Milky Way spiral arm s (e.g.. iLuna et al.l 120061) and 
ring galaxies ( I Vorobvovll2003l) . In regions of low shear, the 
magneto- J eans instability may induce the growth of perturba- 
tions (e.g jKim et~aill2002l) . 



Observational evidence for connection between global disk 
instabilities and s tar formation has been a matter of consider- 
able debate (e.g..lKennicuttlll989HMartin & Kennicuttll2001l: 
Hunter et ail 1 9981; iBoissier et al.ll2003L 120071; iKennicutt et al.l 
20071; lYang et al.ll2007l) . In particular, the existence of global 
star formation thresholds related to Toomre instability (e.g., 
iKennicuttl fl989t iMartin & Kennicuttl 1200 lb has been chal- 
lenged dBoissier et alj |2007Ik It is thus interesting to use 
model galaxies to explore the role of critical surface densi- 
ties in star formation. While proper radiative transfer calcu- 
lations are needed to assign either Ha or UV emissivities to 
the simulated disk galaxies, an indirect comparison between 
simulations and observations can be made by measuring sur- 
face mass densities, velocity dispersions, dynamical proper- 
ties, and a characteristic extent of star formation in the model 
disks. These comparisons are presented in Figure [TTJ where 
the instability parameters £}~ g ' and £ gas /XU are plotted as a 
function of radius in the disk galaxies simulated with the H2D- 
SF+ISRF molecular ISM and star formation model. For com- 
parison with the two-fluid Toomre axisymmetric instability 
criterion, the parameter Q sg is measured by determining the 
locally most unstable wavenumber q from Q s an d Qg with 
a = 0.6 9 to approximate the observatio nal results of Kennicutt 
( 1989) and Ma rtin & Kennicuttl d2001l) . When examining the 
shearing instability criterion, the parameter = 2 is chosen 
for comparison with the observational estimat ion of critical 
surface densities by Mart in & Kennicuttl (1200 ll) and only the 
gas surface density is utilized. A "critical" radius /?sfr con- 
taining 99.9% of the total star formation in each simulated 
disk is measured and used as a proxy for the star formation 
threshold radius observationally determined from, e.g., the 
Ha emission. 

The shearing instability criterion (upper panel Figure [TTJ 
is supercritical (£ gas /£A > 1) over the extent of star forma- 
tion in the low-mass galaxy, but the larger galaxies are either 
marginally supercritical or subcritical. Given the large gas 
content in this dwarf system and that this comparison only 
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Fig. 11. — Stability parameters for two-fluid Toomre instability 
(Gsgi bottom panel) or shear instability (E gas /EA, top panel) as a 
function of radius normalized to the radius Rsfr containing 99.9% 
of the total active disk star formation. Shown are the parameters 
Q7g and T,^ s /T,a for galaxy models with circular velocities of 



= 50km s (squares, light gray lines), 



= 125km s 



(triangles, dark gray lines), and v c i rc = 300 km s~ (circles, black 
lines). The proportionality const ants ap = 0.69 and qa = 2 are 
taken to match the values used by Martin & Kennicutt (2001). For 
each simulated galaxy, the parameter Qsg approaches the transition 
to stability at the radius Rsfr beyond which only very little star 
formation operates. An increase in the stability against growing 
shearing modes also occurs near Rsfr for the simulated disks, but 
the transition is less uniform than for Toomre instability. 



uses the gas surface mass density £ gas , the result is unsur- 
prising. While the v c h- c ~ 125 km s" 1 and v c j lc ~ 300 km s -1 
galaxies are subcritical at R = Rsfr, we note that there is still 
a noticeable drop in £ gas /£A near R = Rsfr ~ 1 for these sys- 
tems. 

Interestingly, with the chosen definition of Rsfr, the bottom 
panel of Figure QT| shows that the two-component Toomre in- 
stability parameter Q sg is a remarkably accurate indicator of 
where the majority of star formation operates in the simulated 
disk galaxies. Each galaxy approaches the Q sg ~ 1 transi- 
tion near Rsfr, demonstrating that much of the star forma- 
tion is operating within Toomre-unstable regions of the disks. 
These results appear most consistent with the obse rvations by 
iKennicuttl d 1 9891) . iMartin & Kennicutil d200 lb . and lYang etal] 
(2007), who find that star formation operates most efficiently 
in regions that are gravitationally unstable. 

Our simulations ar e also qualitativ el y cons istent with the 
theoretical findings of iLi et alJ d2005allbTl2006l) . who used hy- 
drodynamical simulations of isolated galaxies with sink parti- 
cles to represent dense molecular gas and star clusters. They 
found that if the accretion onto the sink particles resulted 
in star formation with some efficiency, the two-component 
Toomre parameter was Q sg < 1 .6 in regions of active star for- 
mation. Although a direct comparison is difficult given their 
use of sink particles, we similarly find that the majority of 
star formation in our simulations occurs in regions of the disk 
where Q sg < 1 . 

In related theoretical work, a gravither mal ins t ability crite- 
rion for star formation was explored by Schaye (2004;), who 
suggested that the decrease in thermal velocity dispersion as- 
sociated with the atomic-to-molecular transition in the dense 
ISM may induce gravitational instability and the correspond- 
ing connection between critical surface densities and star for- 
mation. The simulated galaxies do experience a decline in 
thermal sound speed near the Q sg ~ 1 transition , consistent 
with the mechanism advocated by Schaye (2004). However, 
of our H2D-SF+ISRF simulations (that include physics nec- 
essary for such instability), only the dwarf galaxy model dis- 
plays a sharp transition in Q sg near the critical radius while the 
stability of the other, more massive galaxies are more strongly 
influenced by the turbulent gas velocity dispersion and have 
more modest increases in Q ss ,. Furth er, we note that the 
gravithermal mechanism suggested by Schaye ( 2004) should 
not operate i n isothermal treatments o f the ISM such as those 
presented bv lLi et alJ d2005albl [2006b . Similarly our GD-SF 
simulations, which are effectively isothermal (T w 10 4 K) at 
ISM densities, display similar transitions to 2sg < 1 at the 
critical radius even as such simulations do not include a cold 
phase to drive gravitational instability. 

We stress that the definition of Rsfr may affect conclusions 
about the importance of threshold surface densities. For in- 
stance, if a 90% integrated star formation threshold was used 
Rsfr would be lowered to 0.3-0.6 of the chosen value and 
the subsequent interpretations about the applicability of the 
instability parameter Q s „ could be quite different. Detailed 
radiative transfer calculations may therefore be necessary to 
determine the integrated star formation fraction that brings 
Rsfr in the best agreement with observationally-determined 
critical radii. 



3.6. Temperature and Density Structure of the ISM 

The temperature and density structure of interstellar gas 
is determined by a variety of processes including radia- 
tive heating, cooling, and, importantly, turbulence (see, e.g., 
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FIG. 12. — Temperature structure of the interstellar medium in galaxies of different mass. Shown are the probability distribution functions (PDF) of the 
fractional volume of the ISM with a temperature T (black histograms), for systems with circular velocities of v c j rc = 300km s (left panel), v c j lc = 125 km s 
(middle panel), and v c j rc = 50km s (right panel). The simulations use the FbD-SF+ISRF model to evolve the systems in isolation for t = 0.3Gyr. The peaks in 
the temperature PDF correspond to regions of density-temperature phase space where the cooling time ? coo i has local maxima (i.e., where t coo \ is locally a weak 
function of density). These peaks in the temperature PDF correspond to features in the density PD F, an d the colored arrows highlight regions that can be modeled 
approximately as separate isothermal phases (i.e., lognormals) in the density PDF (see Figure [TSl and j|3.6l of the text for details). 
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FIG. 13. — Density structure of the interstellar medium in galaxies of different mass. Shown are the probability distribution functions (PDF) of the fractional 
volume of the ISM with a density hh (black histograms), for systems with circular velocities of v c j lc = 300km s (left panel), v c j rc = 125 km s (middle panel), 
and v c j rc = 50km s (right panel) for the FbD-SF+ISRF model at t = 0.3Gyr. The features in the density PDF correspond to regions of density-temperature 
phase space where the cooling time ; coo i has local maxima (i.e., where f coo i is locally a weak function of density). The density PDF of the ISM in these simulated 
galaxies is not well-modeled by a single log-normal, except at the high-density tail. The figure shows that the density PDF may be modeled approximately as a 
sum of isothermal phases (i.e., lognormal PDFs in density, colored dashed lines), where the cha racte ristic density of each phase is determined by using the ISM 
equation-of-state to find the typical density of peaks in the temperature PDF (see Figure fTSl and j|3.6l of the text for details). 



Elme greenl 120021) . The calculation of the density proba- 
bility distribution function (PDF) of an isothermal gas, de- 
fined as a volume fraction of gas with a given density, 
suggests that the PDF will resemble a lognormal distribu- 
tion with a dispersion a that scales linearly with the mean 
gas Mach number a oc (3Ai. T he origin of t he l ognor- 
mal form is discussed in detail in iPadoan et al.l d 19971) and 
iPassot & Vazquez-Semadeni ( 1998})- For a non-isothermal 
gas the density PDF may more generally follow a power-law 
distrib ution (IPassot & Vazquez-Semadeni 1998; IScalo et aTl 
[19981) . 

Although the lognormal PD F was initially di s cussed in the 
context of molecular clouds, Wada & Norman (2001) found 
that two-dimensional numerical simulations of a multiphase 
ISM on galactic scales produced a lognormal density PDF 
for dense gas, while at low densities gas followed a flatter 
distribution in log p (see Wad a & Normanll2007l; IWadal 120071 
for more recent investigations of this issue). A similar PDF 



was found i n a cosmq l ogical simulation of galaxy formation 
at z ~ 4 by | Kravtsovl (120031). Motivated in part by the re- 
sult of IWada & Normanl d200ll) , lElmegreenl (T2002) suggested 
that the SK relation could be explained in terms of the density 
PDF of ISM gas if a fraction f c of gas that resides above some 
critical density threshold p c converts into stars with a constant 
efficiency and p c sc ales linea r ly wi th the average ISM den- 
sity (p) . Similarly, iRravtsovl (120031) found that the fraction 
of dense gas in cosmological simulations of disk galaxy for- 
mation at redshifts z > 4 varied with the surface density in a 
manner to produce the Ssfr °c relation at high surface 
densities. The gas density distribution in these simulations 
also followed a log-normal PDF a t large densities and fla t- 
te ned at lower densiti e s (see also lKravtsov & Gne din 2005). 

iJoung & Mac Low! (120061) used three-dimensional models 
of the stratified ISM to demonstrate that supernovae explo- 
sions can act to regulate star formation by inputting energy on 
small scales. Their scheme produced temperature and density 
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structures consistent with those observed in the real ISM, with 
density and temperatu re PDFs that contain multiple peaks. 
Taske r & Brvanl d2007l) found the density PDF in their mul- 
tiphase simulations could be approximated by a single log- 
normal at high densities, but exhibited deviations from the 
lognormal in the form of peaks and dips at smaller densities. 

FigureQ~2]shows the probability distribution function (PDF) 
of the fractional volume in the ISM at a temperature T in 
the H 2 D-SF+ISRF simulations at t = 0.3 Gyr. The temper- 
ature PDF is measured by assigning an effective volume to 
each particle based on its mass and density, counting the par- 
ticles in bins of fixed width in log 10 7\ and then normaliz- 
ing by the total effective volume of the ISM. The tempera- 
ture PDF displays structure corresponding to gas at several 
phases. These temperature peaks correspond to features in 
the density-temperature phase diagram of the ISM, set by the 
cooling and heating processes. The abundance of gas in dif- 
ferent regions of this diagram depend on the local cooling 
time, with regions corresponding to long local cooling times 
becoming more populated. The overall gravitational poten- 
tial of the system also influences the density of gas in the 
disks and thus indirectly influences cooling and abundance 
of cold, dense gas, as denser ISM gas can reach lower equi- 
librium temperatures. Figure [12] shows that the temperature 
PDFs of the galaxies do differ, with more massive galaxies 
containing more volume at low ISM temperatures. In smaller 
mass galaxies, the ISM densities are lower, gas distributions 
are more extended, and the abundance of cold gas is greatly 
decreased. 

Figure [13] shows the PDF of gas densities for each model 
galaxy in the H 2 D-SF+ISRF model at t = 0.3Gyr. The density 
PDFs are measured in the same manner as the temperature 
PDFs, but the particles are binned according to the log 10 «H- 
Each galaxy has a PDF with noticeable features that change 
depending on the galaxy mass scale. None of the density 
PDFs are well represented by a single lognormal distribu- 
tion, reflecting the multiple peaks apparent in the temperature 
PDFs (Figure O. 

Motivated by the multiphase structure of the ISM displayed 
in Figure [T2] a model of the density PDFs can be constructed 
by identifying characteristic temperature peaks in the temper- 
ature PDF and assigning separate lognormal distributions to 
represent each approximately isothermal gas phase. The col- 
ored arrows in Figure [12] identify characteristic temperatures 
in the ISM of each galaxy. The characteristic density corre- 
sponding to each selected temperature, determined by the me- 
dian equation-of-state calculated using Cloudy, is employed 
as the mean of a lognormal distribution used to model a fea- 
ture in the density PDF (colored dashed lines in Figure [T3l. 
The heights and widths of the lognormals are constrained to 
provide the same relative volume fraction in the density PDFs 
as the phase occupies in the temperature PDF, resulting in a 
model shown by the solid gray lines in Figure [T3l 

Clearly the ISM is best understood in terms of a continuous 
distribution of temperatures, represented in the simulations by 
Figure [12] but the multiphase model of the density PDFs pro- 
vides some conceptual insight into the ISM phase structure. 
First, the model of separate isothermal ISM phases can cap- 
ture the prominent features of the density PDFs. Second, the 
galaxies share common features in the density (and temper- 
ature) PDFs, such as a warm neutral phase (log 10 r ~ 3.8) 
and the mix of cold neutral and molecular gas at low tem- 
peratures (log 10 r < 3.5). Third, the relative widths of the 
lognormal distributions representing the dense ISM in each 



galaxy increase with increasing galaxy mass (i.e., the orange 
lognormals centered near log I0 «H ~ 0.5), as do their normal- 
izations. The increasing widths of the lognormals mimic the 
increasing velocity dispersions (and turbulent Mach numbers) 
of the dense gas in the more massive galaxy models. 

Note that for the ISM of the v circ = 300km s" 1 system the 
cold dense phase contributes the majority of the ISM by mass, 
while shock heated gas at T ~ 10 5 K contributes the largest 
ISM fraction by volume. We note that this result qualita- 
tively resembles the effect of the volume -averaged equation- 
of-state utilized in subgrid m odels of the multiphase ISM (e.g. 
ISpringel & Hernquistll2003al) . 

For purposes of calculating the overall SFE of galaxies, our 
results suggest that approximating the density PDF by a single 
lognormal distribution is overly simplistic. The detailed shape 
of the PDF is set by global properties of galaxies and heating 
and cooling processes in their ISM and can therefore be ex- 
pected to vary from system to system. Nevertheless, a single 
PDF modeling may be applicable for the cases when most of 
the ISM mass is in one dominant gas phase (e.g., dense cold 
gas), as is the case for the massive galaxy in our simulations. 
We leave further exploration of this subject for future work. 

4. DISCUSSION 

Results presented in the previous sections indicate that 
star formation prescriptions based on the local abundance of 
molecular hydrogen lead to interesting features of the global 
star formation relations. 

We show that the inclusion of an interstellar radiation field 
is critical to control the amount of diffuse H2 at low gas den- 
sities. For instance, without a dissociating ISRF the low mass 
dwarf galaxy eventually becomes almost fully molecular, in 
stark contrast with observations. We also show that with- 
out the dissociating effect of the ISRF our model galaxies 
produce a much flatter relation between molecular fraction 
fm and pressure P ext , as can be expected from the results of 
Elmeg reenl (Il993ah . Including ^-destruction by an ISRF re- 
sults in a /h2 - Pext relation in excell ent a greement with the ob- 
servations o f IWong & BUtzl (12002b and Blitz & Rosolowskvl 
(f2004l [2006b . 

Our model also predicts that the relation between Ssfr 
and S gas should not be universal and can be considerably 
steeper than the canonical value of « tot = 1 .4, even if the three- 
dimensional Schmidt relation in molecular clouds is univer- 
sal. The slope of the relation is controlled by the dependence 
of molecular fraction (i.e., fraction of star forming gas) on 
the total local gas surface density. This relation is non-trivial 
because the molecular fraction is controlled by pressure and 
ISRF strength, and can thus vary between different regions 
with the same total gas surface density. The relation can also 
be different in the regions where the disk scale-height changes 
rapidly (e.g., in flaring outer regions of disks), as can be see n 
from equation [15] (see also Schaye & Dal la Vecchial 120081) . 
We show that the effect of radial variations in the molecular 
fraction fm and gas scale heights (/isfr and /i gas ) on the SFR 
can be accounted for in terms of a structural, SK-like corre- 
lation, Ssfr oc fmhsFRhZl^gM, that trivially relates the local 
SFR to the consumption of molecular gas with an efficiency 
that scales with the local dynamical time. 

A generic testable prediction of our model is that deviations 
from the SK Ssfr _ S gas relation are expected in galaxies or 
regions of galaxies where the molecular fraction is declining 
or much below unity. As we discussed above, star forma- 
tion in the molecule-poor (fm ~ 0.1) galaxy M33 supports 
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this view as its to tal gas SK power-law index is n tot w 3.3 
dHeveret al.112001 . While our simulations of an M33 ana- 
logue produce a steep SK power-law, the overall SFE in our 
model galaxies lies below that observed for M3 3. However, 
as recently emphasized by iGardan et al.1 (120071) . the highest 
surface density regions of M33 have unusually ef ficient star 
forma tion compared with the normalization of the iKennicutH 
( 1998) relation and so the discrepancy is not surprising. Given 
that dwarf galaxies generally have low surface densities and 
are poor in molecular gas, it will be interesting to examine SK 
relation in other small-mass galaxies. 

Another example of a low-molecular fraction galaxy close 
to home is M31, which has only a fraction fm ~ 0.07 of 
its gas in molecular f orm within 18 kpc of the galactic cen- 
ter ^Nieteneral]|200^). Our model would predict that this 
galaxy should deviate from the total gas SK relation found for 
molecular-rich ga laxies. Observational ly, the SK relation of 
M3 1 measured by Boissi er" ;t al.l (120031) is rather complicated 
and even has an increasing SFR with decreasing gas density 
over parts of the galaxy. Low molecular fractions can also be 
expected in the outskirts of normal galaxies and in the disks 
of low surface brightness gala xies. The latter have molecular 
fracti ons of only f m < 0. 10 (iMatfhews et al.ll2005b iDas et all 
2006), and we therefore predict that they will not follow the 
total gas SK relation obeyed by molecule-rich, higher sur- 
face density galaxies. At the same time, LSBs do lie on the 
same relation between H2 mass and far infr ared luminosity 
as hig her surface brightness (HSB) galaxies (M atthews et al.l 
2005), which suggests that the dependence of star formation 
on molecular gas may be the same in both types of galaxies. 

An alternative formulation of the global star formation re- 
lation is based on the angular frequency of disk rotation: 
Ssfr oc E gas 57. That this relation works in real galaxies is 
not trivial, because star formation and dynamical time-scale 
depend on the local gas density, while ft depends on the to- 
tal mass distribution within a given radius. Although several 
models were proposed to explain s uch a correlation (see, e.g., 
IKennicuttl 1 1 9981: lElm egreen 2002, for reviews), we show in 
33.41 and the second Appendix that the star formation cor- 
relation with il can be understood as a fortuitous correla- 
tion of O with gas density of ft oc p^L., where a « 0.5, for 
self-gravitating exponential disks or exponential disks em- 
bedded in realistic halo potentials. Moreover, we find that 
Ssfr oc Sgasil breaks down at low values of £ gas f2 where 
the molecular fraction declines, similarly to the steepening 
of the SK relation. Our models therefore predict that the 
Esfr oc Eh2^ relation is more robust than the Ssfr oc S gas J7 
relation. 

An important issue related to the global star formation in 
galaxies is the po s sible existence of star form ation thresholds 
(IKennicuttl 1 1998t iMartin & Kennicuttll2001l) . Such thresh- 
olds are expected to exist on theoretical grounds, because the 
formation of dense, star-forming gas is thought to be facil- 
itated by either dynamical insta bilities (see, e.g., lElmegreenl 
2002; McKee & Ostriker 200 71 for compr ehensive reviews) 
or gravithermal instabilities dSchavell2004l) . We find that the 
two-component Toomre instability threshold that accounts for 
both stars and gas, Q sg < 1, works well in predicting the tran- 
sition from atomic gas, inert to star formation, to the regions 
where molecular gas and star formation occur in ou r simu- 
lations. Our re sults are in general agreement with iLi et al.1 
(2005b a, 2006), who used sink-particle simulations of the 
dense ISM in isolated galaxies to study the relation between 
star formation and the development of gravitational instabil- 



ities, and with observations of star formation in the LMC 
(Yang et al. 2007). Our simulations further demonstrate the 
importance of accounting for all mass components in the disk 
to predict correctly which regions galactic disks are gravita- 
tionally unstable. 

Given that our star formation prescription is based on 
molecular hydrogen, the fact that Q sg is a good threshold in- 
dicator may imply that gravitational instabilities strongly in- 
fluence the abundance of dense, molecular gas in the disk. 
Conversely, the gas at radii where the disk is stable remains 
at low density and has a low molecular fraction. We find 
that in our model galaxies, the shear instability criterion of 
Elmegreen ( 1993a) does not work as well as the Toomre Q sg - 
based criterion. Almost all of the star formation in our model 
galaxies occurs at surface de nsities S gas > 3 M Q pc" 2 , which 
is formally consistent with the lSchay e (2004) constant surface 
density criterion for gravithermal instability. However, as Fig- 
ure [6] shows, we do not see a clear indication of threshold at 
a particular surface density and our GD-SF models that have 
an effectively isothermal ISM with T w 10 4 K (and hence do 
not have a gravithermal instability) still show a good corre- 
lation between regions where Q sg < 1 and regions where star 
formation operates. 

Our results have several interesting implications for inter- 
pretation of galaxy observations at different epochs. First, 
low molecular fractions in dwarf galaxies mean that only a 
small fraction of gas is participating in star formation at any 
given time. This connection between SFE and molecular hy- 
drogen abundance may explain why dwarf galax ies are still 
gas ri ch today compared to larger mass galaxies (Ge ha et all 
2006), without relying on mediation of star formation or gas 
blowout by supernovae. Note that a similar reasoning may 
also explain why large LSBs at low redshifts are gas rich but 
anemic in their star formation. Understanding the star forma- 
tion and evolution of dwarf galaxies is critical because they 
serve as the building blocks of larger galaxies at high red- 
shifts. Such small-mass galaxies are also expected to be the 
first objects to form large masses of stars and should there- 
fore play an important role in enrichment of primord ial gas 
and the cosmic star formation rate at high redshifts (iHopkinsI 
2004: lHopkins & Beacomll2006l) . 

The star-forming disks at z ~ 2 that may be progen- 
itors of low-redshift spiral galaxies are o bserved to lack 
centra lly-concentrated bulge components (Elmegre en et al.1 
120071) . Given that galaxies are expected to undergo frequent 
mergers at z > 2, bulges should have formed if a signifi- 
cant fraction of baryons are conv erted into stars during such 
mergers ("e.g.. lGnedin et al.ll2000l) . The absence of the bulge 
may indicate that star formation in the gas rich progenitors 
of these z ~ 2 systems was too slow to convert a signifi- 
cant fraction of gas into stars. This low SFE can be under- 
stood if the high cosmic UV background, low-metallicities, 
and low dust cont ent of high-z gas disk s keep their molec- 
ular fractions low (Pelup essy et al.ll2006l) . thereby inhibiting 
star formation over most of gas mass and keeping the pro- 
genitors of the star-forming z ~ 2 disks mostly gaseous. Gas- 
rich progenitors may also help explain the prevalence of ex- 
tended disks in low-redshift galaxies despite the violent early 
merger histories characteristic of ACDM universes, as gas- 
rich mer gers can help build high-angular momentum disk 
galaxies (Robertson et al. 2006a). Mergers of mostly stellar 
disks, on the other hand, would form spheroidal systems. 

Our results may al so provide insight into the interpreta- 
tion of the results of IWolfe & Che n (2006), who find that 
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the SFR associated with neutral atomic gas in damped Ly- 
man alpha (DLA) syst ems is an order of magnitude lower than 
predicted by the local iKennicuttl JT998) relation. The DLAs 
in their study sample regions with column densities A^h ~ 
2 x 10 21 cm" 2 , or surface gas densities of Eqla ~ 20 M Q pc~ 2 , 
assuming a gas disk with a thickness of h w 100 pc. Suppose 
the local SK relation steepens from the local relation with the 
slope «o = 1 .4 to a steeper slope n \ below some surface den- 
sity Eb > Sdla- Then for E gas < £],, the SFR density will 
be lower than predicted by the local relation by a factor of 
(Egas/Sb)" 1 "" - For n = 1.4 and n x = 3 the SFR will be sup- 
presse d by a factor of > 10 for E gas /Eb < 0.25. Thus, the 
results [Wolfe & Chenl (120061) can be explained if the total gas 
SK relation at z ~ 3 steepens below E gas < 100 M©pc~ 2 . We 
suggest that if the majority of the molecular hydrogen at these 
redshifts resides in rare, compact, and dense systems (e.g., 
IZwaan & Procha ska 2006), then both the lack of star forma- 
tion and the rarity of molecular hydrogen in damped Lya ab- 
sorbers may be explained simultaneously. 

Our results also indicate that the thermodynamics of the 
ISM can leave an important imprint on its density probabil- 
ity distribution. Each thermal phase in our model galaxies 
has its own log-normal density distribution. Our results thus 
imply that using a single lognormal PDF t o build a model 
of global star formatio n in galaxies (e.g.. lElmegreenll2002t 
Wada & Norman 2007) is likely an oversimplification. In- 
stead, the global star formation relation may vary depend- 
ing on the dynamical and thermodynamical properties of the 
ISM. We can thus expect differences in the SFE between the 
low-density and low-metallicity environments of dwarf and 
high-redshift galaxies and the higher-metallicity, denser gas 
of many large nearby spirals. Note that many of the results 
and effects we discuss above may not be reproduced with 
a simple 3D density threshold for star formation, as com- 
monly implemented in galaxy formation simulations. Such 
a threshold can reproduce the atomic-to-molecular transition 
only crudely and would not include effects of the local inter- 
stellar radiation field, metallicity and dust content, etc. 

A clear caveat for our work is that the simulation resolution 
limits the densities we can model correctly. At high densities, 
the gas in our simulations is over-pressurized to avoid numer- 
ical Jeans instability. The equilibrium density and tempera- 
ture structure of the ISM and the molecular fraction are there- 
fore not correct in detail. Note, however, that our pressuriza- 
tion prescription is designed to scale with the resolution, and 
should converge to the "correct" result as the resolution im- 
proves. In any event, the simulations likely do not include all 
the relevant physics shaping density and temperature PDFs of 
the ISM in real galaxies. The results may of course depend on 
other microphysics of the ISM as they influence both the tem- 
perature PDF and the fraction of gas in a high-density, molec- 
ular form. Future simulations of the molecular ISM may need 
to account for new microphysics as they resolve scales where 
such processes become important. 

5. SUMMARY 

Using hydrodynamical simulations of the ISM and star for- 
mation in cosmologically motivated disk galaxies over a range 
of representative masses, we examine the connection between 
molecular hydrogen abundance and destruction, observed star 
formation relations, and the thermodynamical structure of the 
interstellar medium. Our simulations provide a variety of new 
insights into the mass dependence of star formation efficiency 
in galaxies. A summary of our methodology and results fol- 



lows. 

1. A model of heating and cooling processes in the in- 
terstellar medium (ISM), including low-temperature 
coolants, dust heating and cooling processes, and heat- 
ing by the cosmic UV background, cosmic rays, and 
the local interstellar radiation field (ISRF) , is calculated 
using the photoionization code Cloudy dFerland et alJ 
1998). Calculating the molecular fraction of the ISM 
enables us to implement a prescription for the star for- 
mation rate (SFR) that ties the SFR directly to the 
molecular gas density. The ISM and star formation 
model is implemented in the SPH/N-body code GAD- 
GET2 (ISpringelll2005l) and used to simulate the evolu- 
tion of isolated disk galaxies. 

2. We study the correlations between gas surface den- 
sity (Egas), molecular gas surface density (Eh2)> and 
SFR surface density (Esfr). We find that in our most 
realistic model that includes heating and destruction 
of H2 by the interstellar radiation field, the power 
law index of the SK relation, Esfr oc E"™, (measured 
in annuli) varies from n tot ~ 2 in massive galaxies 
to n to t > 4 in small mass dwarfs. The correspond- 
ing slope of the Esfr oc E^' 1 molecular-gas Schmidt- 
Kennicutt relation is approximately the same for all 
galaxies, with n mo \ w 1.3. These results are consis- 
tent with obse rvations of star f ormation in different 
galaxies (e.g. , IKennicuttl 119981: IWong & Blitzl 120021: 
Boissier et all 120031: [Hever et all l2004t iBoissier etail 
2007tlKennicutt et alJ]2007l) . 

3. In our models, the SFR density scales as Esfr 
/H2/isFR/ig;l s 5 Ej (a 5 s , where h m is the scale-height of the 
ISM and /?sfr is the scale-height of star-forming gas. 
The different Esfr _ E gas relations in galaxies of differ- 
ent mass and in regions of different surface density in 
our models therefore owe to the dependence of molec- 
ular fraction /h2 and scale height of gas on the gas sur- 
face density. 

4. We show that the Esfr oc E gas il and Esfr oc Eh2^ cor- 
relations describe the simulations results well where 
the molecular gas and total gas densities are compara- 
ble, while the simula tions deviate from Esfr oc E gas i7 
("e . g . . lKennicutdll998l) at low E gas owing to a declining 
molecular fraction. We demonstrate that these relations 
may owe to the fact that the angular frequency and the 
disk-plane gas density are generally related as il oc y/p 
for exponential disks if the potential is do minated by 
either the disk, a iNavarro et al.l Jl996) halo, Hernquist 
( 1990) halo, or an isothermal sphere. The correlation of 
Esfr with O is thus a secondary correlation in the sense 
that oc J~p is set during galaxy formation and £1 does 
not directly influence star formation. 

5. The role of criti cal sur f ace de nsities for shear instabil- 
ities (E^) and|Toomre (1964) instabilities (Eg) in star 
formation (e.g. . [Martin & Kennicuttll2~001h is examined 
in the context of the presented simulations. We find that 
the two-component Toomre instability criterion Q sg < 1 
is an accurate indicator of the star-forming regions of 
disks, and that gravitational instability and star forma- 
tion are closely related in our simulations. Further, the 
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Q sg criterion works even in simulations in which cool- 
ing is restricted to T > 10 4 K where gravithermal insta- 
bility cannot operate. 

6. Our simulations that include H2-destruction by an ISRF 
naturally reproduce the observed scaling fm oc P^ be- 
tween molecular fra ction and external p ressur e (e.g ., 
IWong & Blit jl200l iBlitz & Rosolowskvll200l [2006), 
but we find that simulations without an ISRF have a 
weaker scaling f H 2 oc P^ 4 . We calculate how the con- 
nection between the scalings of the gas surface density, 
the stellar surface density, and the ISRF strength influ- 
ence the /h2 - Pext relation in the ISM, and show how 
the simulated scalings reproduce the fm~Pext relation 
even as the power-law index of the total gas Schmidt- 
Kennicutt relation varies dramatically from galaxy to 
galaxy. 

7. We present a method for mitigating numerical Jeans 
fragmentation in Smoothed Particle Hydrodynamics 
simulations that uses a density-dependent pressuriza- 
tion of gas on small scales to ensure that the Jeans 
mass is properly resolved, simila r to techniques used 
in grid-based simula tions (e.g., iTruelove et al.l 119971: 
Machacek et al. 2001). The gas internal energy u at the 

-2/3 

Jeans scale is scaled as u oc 'w Jeans , where mj eans is the lo- 
cal Jeans mass, to ensure the Jeans mass is resolved by 
some Njeans number of SPH kernel masses 2Af ne i g h'wsPH, 
where Af ne i g h is the number of SPH neighbor particles 
and mspH is the gas particl e mass. For the simul ations 
presented here, we find the Bate & Burkerll (1 19971) cri- 
terion of Afjeans = 1 to be insufficient to avoid numer- 
ical fragmentation and that Nj eans ~ 15 provides suf- 
ficient stability against numerical fragmentation over 
the time evolution of our simulations. Other simula- 
tions may have more stringent resolution requirements 
(e.g.. ICommercon et al.l 12008b . We also demonstrate 
that isothermal galactic disks with temperatures of T = 
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10 4 K may be susceptible to numerical Jeans instabili- 
ties at resolutions common in cosmological simulations 
of disk galaxy formation, and connect this numerical ef- 
fect to possible angular momentum deficiencies in cos- 
mologically simulated disk galaxies. 

The results of our study indicate that star formation may 
deviate significantly from the relations commonly assumed 
in models of galaxy formation in some regimes and that these 
deviations can be important for the overall galaxy evolution. 
Our findings provide strong motivation for exploring the 
consequences of such deviations and for developing further 
improvements in the treatment of star formation in galaxy 
formation simulations. 
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APPENDIX 



iTruelove et al.l (1 19971) showed that hydrodynamical calculations that do not properly resolve pressure forces on scales where 
self-gravity begins to dominate over gaseous p r essure support are susceptible to a numerical instability that leads to the rapid 
artificial fragmentation of gas. Truel oveet al.l ([1997) suggest a constraint on the spatial resolution of hydrodynamical mesh 
calculations such that the Jeans length Aj eans is locally larger than the grid resolution Ax. This resolution requirement for grid 
codes can be expressed in terms of a constraint on the "Jeans number" J = Ax/Xj ew s < 0.25 to properly resolve the Jeans scale 
and prevent such artificial fragmentation. 

Similarly, the inability of SPH codes to capture pressure forces properly on scales less than the smoothing length can also allow 
fo r numerically-induced f ragmentation if the smoothing length is of order the Jeans length dBate & Burkert fl997l) . As discussed 
by Bate & Burkert ( 1997) and in 32.31 the requirement that the Jeans scale is resolved an SPH calculation by some number A?j eans 
of SPH kernels is 27Yj eans 7Y ne ighm gas < ra Jeans , where m gas is the SPH particle mass and jV n ri g h is the number of SPH neighbors. 
A similar requirement can be formulated in terms of relevant spatial scales (e.g.. iKlein et al.ll2004l) . Owing to the temperature 
dependence of the sound speed (c, oc \ff for ideal gas) and the dependence of the Jeans scale on density and temperature, it is 
clear that for a given number of gas particles 7Y gas used in an SPH calculation, the lower the temperature floor allowed the more 
potentially susceptible the simulated fluid will be to numerical Jeans fragmentation. 

Simulations of the ISM in galaxies with molecular coolants must be especially cautious given the typical numerical resolutions, 
high densities, and potentially low temperatures. For example, for the typical conditions of the average ISM in a Milky Way- 
size galaxies (number density of «h = 1cm" 3 , sound speed c s = lOkms" 1 ), the Jeans mass is mj ean s ~ 6.5 x 10 7 M Q , which 
is comparable to or smaller than the gas mass resolution in many recent simulations of disk galaxy formation. However, for 
molecular gas («h ^ 100cm" 3 ) with a sound speed c s = lkms" 1 (corresponding to the temperature T = 100 K), the Jeans mass 
would be only mj ean s ^ 6.5 x 1O 3 M0 — far smaller than the gas mass resolution of the largest MW-sized galaxy formation 
simulations performed to date. 

For SPH calculations in the literature, typical approaches for addressing numerical Je ans instability include improving the 
numerical resolution such that the Jeans scale is resolved with Afj eans = 1, as suggested bvlBate & Burkertl dl997l) . or using sink 
particles to represent high densi ty regions of the ISM (e.g. iLi et al.i r2005b, 2006; Booth et al. 2007). Each of these approaches 
have possible drawbacks, as the iBate & Burkertl (1 19971) resolution criterion can clearly depend on the density, depending on the 
gas equation-of-state, and the sink particle method, depending on its implementation, may neither avoid numerical fragmentation 
nor properly capture the com plicated chemical, thermodynamical, and radiative processes that operate in dense ISM gas (see, 
e.gJCommercon eta l. 2008). Another possible solution is to use a stiffer equation-of-state for the gas (e.g. JSpringel et aLll2005l : 
ISchave & Dalla Vecchiall2008l) . 

Instead, the simulations pres ented in this paper ado pt a pressure floor on the scale of the Jeans length through an approach 
similar to that implemented bv lMachacek et all (|2001) in grid codes. iMachacek et al.l yOOl) adopt an effective pressure P e ff on 
the scale of the finest grid resolution Ax, given by 



where pb is the baryon density in the cell, p is the mean molecular mass, G is the gravitational constant, and K sa 100 is a 
parameter that sets the effective resolution of their simulations. If the effective pressure P oc pu on the scale of the SPH kernel 
is set to ensure the Jeans scale is over-resolved by a factor Afj ean s, then the internal energy of gas particles with kernel mass 
mspH ~ ^V n eigh'Wg a s should be scaled as 



which is Equation [Tol from {32.31 A factor 7Yj eans = 1 would correspond to a d ensity-dependent version of the Bat e & Burkertl 
d 19971) criterion , while applying the effective pressure of Macha cek et al.l d2001l) with K ~ 100 and equating Ax ~ h would imply 
Njeans > 1000. IKlein et all d2004l) provide a survey of the literature and find that values of Nj ems ~ 0.4-30 for a fixed maximum 
density are in common usa ge. 

As briefly mentioned in 32.31 we have performed a variety of resolution studies where we examine the range of iVj eans = 1-100 
and found A^j ea ns ~ 15 to provide sufficient stability for the time interval in our simulations. Note, however, that more stringent 
criteria are likely to be required for longer periods of evolution for galaxies in cosmological simulations. The resolution studies 
utilized a galaxy model initialized with a virial velocity V v ; r = 160km s" 1 , spin A = 0.033, NFW concentration c = 9, disk mass 
fraction ma = 0.05, bulge mass fraction m\, = 0.0167, and disk gas fraction / gas = 0.2 according to the method described in 32.61 
The dark matter halo, bulge, stellar disk, and gaseous disk were resolved with Ndm = 120,000, N\, = 40,000, Nd = 128,000, and 
N s = 400,000 particles respectively, and the halo and baryonic particle gravitational softenings were set to 6dm = 100/r'pc and 
eb = 50/r'pc. The mass contained within the SPH kernel with A^eigh = 64 for this galaxy model was typically wsph = 2 x 10 6 M Q . 
The galaxy model was evolved with two ISM models, the atomic cooling treatment with a temperature floor of T = 10 4 K (the 
GD-SF model) and the atomic and molecular cooling treatment with a temperature floor of T = 10 2 K (the H2D-SF model). 
Figure [T4l shows the disk galaxy evolved with a density-dependent pressure floor described by Equation [TOl with Nj ams =15 
(upper row) and without a pressure floor (lower row) evolved in isolation for t = 0.2Gyr. The systems without a pressure floor 
quickly fragment if the effective temperature floor is either T = 10 2 K (lower middle panel) or T = 10 4 K (lower left panel), as the 
Jeans mass is initially not resolved at the highest gas densities for either temperature. By design, the simulations that include the 
suggested pressure floor do not fragment on the Jeans scale (upper row). The pressurization on the Jeans scale used to ameliorate 



P et f = KGp 2 b Ax 2 /p 



(1) 




(2) 
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Fig. 14. — Results of simulations of an isolated galaxy model with a circular velocity of v c j rc = 160km s _1 evolved with (upper row) and 
without (lower row) a pressure floor to resolve the local Jeans mass and length. The image intensity tracks the surface density while the color 
reflects the gas effective temperature. For galaxy systems of this mass and disk structure with gas fractions of / gas = 0.2 and simulated with 
N = 400, 000 gas particles, the local Jeans mass in dense regions is smaller than the typical mass in an SPH kernel if ISM is modeled with 
either atomic cooling and a T = 10 4 K temperature floor (lower left panel) or atomic and molecular cooling with a T = 10 2 K temperature floor 
(lower middle panel). These systems are susceptible to numerical fragmentation since, owing to the limi ted resolution, the lo cal gas pressure 
is artificially low on the scale where the gas is marginally stable against collapse through self-gravity dTruelove et aljfl997n . To avoid this 
numerical artificiality, a pressure-floor is introduced to increase the local sound speed until the Jeans scale is resolved by several SPH kernels. 
Our prescription allows for galaxy models with either atomic cooling at temperatures T > 10 4 K (upper left panel) or atomic and molecular 
cooling at temperatures T > 10 2 K (upper middle panel) to be evolved stably. Systems that are extremely gas-rich (e.g., / gas = 0.99, right panels) 
can still become globally unstable and undergo large scale fragmentation if a pressure floor is used (upper right panel), but fragmentation occurs 
at much smaller scales and is qualitatively different without such preventative measures (lower right panel). 

artificial fragmentation does not insure that simulated disks will avoid other instabilities, such as global Toomre ( 1964) instability. 
The right panels of Figure H~4l show the W = 160km s" 1 disk model with the disk gas fraction increased to / gas = 0.99, evolved 
for t = 0.2Gyr in isolation. The gas rich model simulated without a pressure floor violently fragments, but the gas rich model 
evolved using a pressure floor a lso u ndergoes fragmentatio n on a much larger mass scale. As discussed in Springel & Hernquist 
(2003a), Robert son et al.l (120041) and Springel et al. (2005), systems like this gas rich model that are violently Toomre-unstable 
when using an ideal gas or weakly pressurized equation-of-state can be stabilized by adopting a more strongly pressurized ISM 
equation-of-state, but we do not examine such models here. 

Separating artificial fragmentation and physical Jeans fragmentation owing to self-gravity may be difficult without detailed 
resolution studies that examine the convergence of the growth rate of perturbations on different size scales. For Galaxy-mass 
systems with the disk structure we utilize and gas fractions of / eas = 0.2 simulated with N = 400,000 gas particles, the local Jeans 
mass in dense regions is smaller than the typical mass in an SPH kernel if either an ISM model with atomic cooling and a T = 10 4 K 
temperature floor or (the GD-SF model) or an ISM model with atomic + molecular cooling and a T = 10 2 K temperature floor (the 
H2D-SF model) is used, by a factor of »Jj eans /2A'jeans'WsPH ^100 for a T = 10 4 K temperature floor or mj ea ns/2A(ieans'«sPH > 10,000 
for a T = 10 2 K temperature floor. The simulations that include a pressure floor maintain TOj ea ns/2./Vj ea n S ?wspH ~ 1 by design, 
with a small scatter induced by fluctuations in wsph from changing particle masses from star formation and the small allowed 
range of neighbors A^ e i g h = 62 - 66. Given the values of the Jeans mass for the typical densities and temperatures of the ISM 
('Mjeans ~ 10 3 - 10 7 M Q , see above), our calculations thus suggest that simulating molecular cooling in Milky Way-sized disks of 
M gas as 10 10 M© with SPH without using an artificial pressure floor would require N ~ 10 9 - 10 10 gas particles to avoid numerical 
Jeans fragmentation if gas is allowed to cool to 100 K. We delay such resolution studies for future work (and faster computers). 
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Artificial fragmentation of the ISM in disks may be cl osely related to the angu lar momentum problem for producing realistic 
disk galaxies in cosmological simulations (see, e.g., §6 of Robertson et al.l (120041) and Robertson et al. (2006a)). The simulations 
of a Milky Way analogue galaxy model presented in this appendix contain N s - d s = 400,000, which ex cee ds the numerical r e solu- 
tion of many disk galaxies formed in cosmological simulations to date. Both Truelove et all (1 19971) and lBate & Burkertl (1 19971) 
comment on the possible implications of numerical Jeans fragmentation for studies of galaxy formation; based on these test simu- 
lations we suggest that to reliably follow the formation and evolution of disks using SPH simulations a pressure floor similar to our 
Equation [TOl should be used. Prescriptions th at mimic pressurization from the multiphase structure of the ISM that may develop 
through feedback from star formation (e.g..|Yepes et alJI 19971 iHultman & Pharasvn 1999; Springel 2000; Thac ker~& Couchmanl 
2000; Springel & Hernquist 2003a; IStinson et alj|2006t ISchave & Dalla Vecch ia 2008) may also serve to ameliorate numerical 
Jeans fragmentation by raising the sound sp eed as a function of t he mean ISM density. This affect is apparent in numerical tests 
examining the stability of isolated disks by Springel et al. (2005) using various ISM equations-of-state. We plan to revisit the 
implications of numerical Jeans fragmentation on cosmological simulations of gala xy formation in future work . Lastly, we note 
that since disks that are Toomre-unstable are also roughly Jeans-unstable (see, e.g.. lKrumholz & M cKee 2005), numerical Jeans 
fragmentation will artificially reinforce any natural correspondence between gravitational instability and star formation. Sim- 
ulations that examine the relation between star formation and the Toomre instability parameter Q sg should therefore explicitly 
demonstrate that the Jeans scale is properly resolved at all densities. 
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As discussed in 93.41 the simulated galaxies display a correlation between the SFR surface density Esfr and the product of 
the total gas surface mass density S gas and the angular frequency ft in regions where the molecular gas and t otal gas densities 
are comparable. The correlation in our simulations is similar to the observed £sFR-£ gas f2 correlation (e.g., Ken nicuttlll998t 
iKennicutt et alJl2007h . Since the star formation relation implemented in the simulations calculates the SFR density p\ as pro- 
portional to the product of the molecular fraction fui and the gas density scaled with the local dynamical time p 1 5 (Equation 
@), the Esfr - S gas il correlation would occur if the angular frequency scaled with the density as il oc Given that the gas 
density is a local quantity while orbital frequency is determined by the overall mass distribution and potential of the galaxy, it is 
not immediately clear why such a relation should hold. Here we show that the relation between p g and with a similar scaling 
exists for self-gravitating exponential disks or exponential disks embedded in realistic dark matter potentials. 

Given the axisymmetry of the disk galaxy models, the angular frequency is related to the density through Poisson's equation 

d 2 $ d 2 $ , 

d? + dF + " 4 ^ GAotal - (1) 

(Equation[26]in 93.4b . Below, the scaling of the angular frequency ft with the local disk density p^ is calculated in the limiting 
cases of potentials dominated by disks (p t otai ~ Pdisk) or dark matter halos (p t otai ~ Pdm) for exponential disks (paisk oc exp[-R / /?</]), 
and we show that for the examined cases f2 oc ^Pdi?kB{R/Rd), where B(R/Rd) ~ 0(1) is a weak function of the disk radius R that 
depends on the form of the dominate potential. 

Exponential Disk-Dominated Potentials, p^ » /?dm 

For a flattened system, iBinnev & Tremainel d!987l) recognized that as the vertical scale h of the density distribution decreases 
the first term on the left hand side of Equation[26]begins to dominate and eventually 

d 2 $> 

-^^Gp. (2) 
If the surface density of the disk is exponential, the axisymmetric density distribution of the disk can be written as 

A£n 

Pd(R, z)=~^ exp(-R/R d \f(z) (3) 
Rd 

where f(z) oc sech 2 (z) is determined by Equation[2] and the constantA is determ ined by taking /(0) = 1 and the integral constraint 
from the disk mass. Using the methods of §2.6.3 of IBinnev & Tremainel (11987). one can show that the rotation curve in the plane 
of a very thin exponential disk is 

v 2 c (R) = 4nGZ R d y 2 [7o(y)J5>(y)-/i(y)*iO0] , (4) 

where y = R/Rd dFreemanll 19701) . Substituting for T,o/Rd from Equation[3] as is appropriate if the gas densities in the disk plane 
decrease in proportion to the surface density, we find 

Odi sk (l?)= {47rGA- I p d (/?,z)exp( 3 ;) x [7 o (y)^()0-/i(y)*iO0j} 1/2 (5) 
= sfpd~B^(R/R d ). (6) 

The function BdiskO*) is a weak function of the disk surface density £(/?). Over the inner 95% of the disk mass (x ~ 
{0.355,4.75}), the fractional deviation from the fi a ^fpd correlation, normalized to the disk properties at one scale length, 
is only B(0.355)/B(l)- 1 = 0.48 near the disk center and B(4.75)/B(l)- 1 = -0.26 in the disk exterior. 
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Halo-Dominated Potentials, poM >> Pdisk 
In the case of jO to tai ~ Pdm, the angular frequency can be determined directly from the integrated mass M within a radius R as 

GM« R) 



51- 



Ri 



(7) 



There are a wide variety of possible dark matter halo profiles, but we limit ourselves to t hree s alient models. 

Dark matter halos in cosmological simulations on average follow the Nava rro et al.l (jl996) (NFW) dark matter halo profile 
given by 



(VR T +?/r s )(l + VR T +?/r s ) 2 
where p s is the density at the scale radius r s . If we substitute for the disk density 

p s = fl x p d (R,z)exv(R/R d -r s /R d ), 
then the angular frequency in the disk plane of a NFW profile can be written 

^NFW(^) = VT^/^NFW [R/ r s ) , 

where the function 



B NFW (x) = < \6irGf d exp 



r s (x- 1) 
Rd 



(l+x)ln(l+x)-jc 



x 3 (l+x) 



1/2 



(8) 

(9) 
(10) 

(11) 



is, again, a weak function of the disk surface density £(/?). Over the inner 95% of the disk mass, the fractional deviation from 
the 51 oc ^fp~2 correlation, normalized to the disk properties at one scale length, is only Z?NFw(0.355^/r s )/BNFw(^/'"j) _ 1 = 0.27 
near the disk center and B^p\sf(4.75 R d /r s ) /B^pw(R d jj s) - 1 = 1 .45 in the disk exterior. 
Alternatively, we can consider the IHernquistl (I 1 9901) profile 



(12) 



( VW+e/a){ 1 + VW+e/af ' 

where p a is the density at a scale radius a. The IHernquistl (1 19901) profile is used as the DM halo profile in the galaxy models 
presented in this paper and has a density and potential similar to the NFW profile within the scale radius. Substituting for the 
disk density, we can write 

p a = fd 1 Pd(R,z)exp(R/R d -a/R d ). (13) 
The mass profile of the IHernquistl (1 1 990 ) halo provides 

n H (R)=^plB H (R/a) (14) 

where the function 



B H (x)= 32 7 rG/7 1 exp 



a(x— 1) 



Rd 



1 



x(l+x) 2 



1/2 



(15) 



has a scaling nearly identical to the function BnfwM when x << 1, and is a correspondingly weak function of the disk surface 
density £(/?). Over the inner 95% of the disk mass, the fractional deviation from the 51 oc ^fpd correlation, normalized to the disk 
properties at one scale length, is only Bu{0.355R d / a) / Bn{R d / a)- 1 = 0.26 near the disk center and B\i(A.15Rd / a) / Bn(R d / a)- \ = 
1 .48 in the disk exterior. 

In real gala xies, it is thought that th e dark matter halo will adiabatically compress in response to the dissipational formation of 
the disk (e.g., Blumenthal et al. 1986). Models of adiabatic compression s uggest that the dar k matter halo response tends to drive 
the halo profile to resemble an isothermal sphere over much of the disk (Gne din et alj |2()04). The isothermal density profile can 
be written 

ID \ P° 

p iso (R,z)^ 



(16) 

(17) 

(18) 
(19) 

The function Bi so (x) is also weak function of the disk surface density £(/?), and is of the same order as Z?nfw or Z?h- Over the 
inner 95% of the disk mass, the fractional deviation from the 51 oc ^fpd correlation, normalized to the disk properties at one scale 
length, is only B iso (0.355)/B iso (l)- 1 = 1.04 near the disk center and B iso (4.75)/B iso (l)- 1 = 0.37 in the disk exterior. 



(VR T Tz 1 /r ) 2 ' 

where po is the density at an arbitrary scale radius r . For simplicity, we set ro = Rd and substitute for the disk density as 

Pa = r d l pd(R,0)^(R/R d -\). 
Integrating the isothermal density profile to find the enclosed mass provides the angular frequency in the disk plane as 

O iso (tf) = ^Pd~B iso (R/R d ) 

where 



Bi so (x) = {4^G/7V 2 exp(x- 1)} 



1/2 
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FIG. 15. — Relation between surface density S and angular frequency Q, for a variety of disk galaxy mode l potentials. Shown is the qua ntity S 1 5 tx Ssfr 
plotted against the product S£7 for exponential disks in a disk-dominated (solid line), Navarro et al. 1 1997) halo-dominated (dotted line), Hernquist 1 1990) 
halo-dominated (dashed line), or isothermal halo-dominated (dash-dotted line) potential, normalized to the disk properties at an exponential scale length. 
Deviations from the S 1 * 5 ex T,Q correlation (indicated with a dashed grey line) are typically ~ 0(1) and are largest in the very center and exterior of the disk. 
The calculation suggests that the E$pr oc T,Q may reflect the scaling of the local density in exponential gas disks with the characteristic shape of the potentials 
that host them, rather than an a fundamental connection between star formation and global disk processes. 



The Ssfr — Correlation 

Figu re [151 summarizes the behavior of the functions i7disk, ^nfw, ^h, and f2j so . In the iHernquistl ( 1 19901) and iNavarro et al.l 
( 1996) halo potential models, the disk scale length-to-halo scale radius ratios were set to Rd/r s = 0.1 appropriate for observed 
spiral galaxies. The quantity £(/?)' 5 oc Ssfr is plotted as a function of the product Y,(R)tt(R), normalized to the disk properties 
at the scale length Rd over the radial range R = [0,5Rd], for the potentials calculated in this Appendix. The quantity ££1 closely 
tracks E 1-5 throughout the disk, demonstrating that deviations from the il cx ^fpd correlation are typically ~ O(l). Hence, the 
Esfr _ £f2 relation may be understood in terms of the scaling of the local density in exponential disks with the angular frequency 
of realistic galaxy potentials, rather than some more fundamental relation connecting star formation processes in the ISM with 
global disk properties. Regions of galaxies with large molecular gas fractions (Eh2 ~ S gas ) are then expected to follow the 
EsFR-£f2 relation, as demonstrated by our simulations (see 33.41) . We note that while a different normalization would reduce 
the typical deviation from the £ scaling, the high-density portion of the relation corresponding to the region where the rotation 
curve is rapidly increasing should produce a s hallower relation than Ssfr oc Xfl Weak evidence for this may exist in the available 
data (e.g., the inner-disk region of Figure 9 of Ken nicutt et ai1l2007l) . 



